Skip to content

Flu Core API Reference

Docstrings and references for flu_core module.

base_path = clt.utils.PROJECT_ROOT / 'flu_instances' / 'texas_input_files' module-attribute

AbsoluteHumidity

Bases: Schedule

Source code in CLT_BaseModel/flu_core/flu_components.py
class AbsoluteHumidity(clt.Schedule):

    def __init__(self,
                 init_val: Optional[np.ndarray | float] = None,
                 timeseries_df: pd.DataFrame = None):
        """
        Args:
            init_val (Optional[np.ndarray | float]):
                starting value(s) at the beginning of the simulation
            timeseries_df (Optional[pd.DataFrame] = None):
                must have columns "date" and "absolute_humidity" --
                "date" entries must correspond to consecutive calendar days
                and must either be strings with `"YYYY-MM-DD"` format or
                `datetime.date` objects -- "value" entries correspond to
                absolute humidity on those days. Identical to
                `FluSubpopSchedules` field of same name.
        """

        super().__init__(init_val)

        self.timeseries_df = timeseries_df

    def update_current_val(self, params, current_date: datetime.date) -> None:
        self.current_val = self.timeseries_df.loc[current_date, "absolute_humidity"]

    def postprocess_data_input(self) -> None:
        self.timeseries_df = self.timeseries_df.set_index('date')

__init__(init_val: Optional[np.ndarray | float] = None, timeseries_df: pd.DataFrame = None)

Parameters:

Name Type Description Default
init_val Optional[ndarray | float]

starting value(s) at the beginning of the simulation

None
timeseries_df Optional[pd.DataFrame] = None

must have columns "date" and "absolute_humidity" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date objects -- "value" entries correspond to absolute humidity on those days. Identical to FluSubpopSchedules field of same name.

None
Source code in CLT_BaseModel/flu_core/flu_components.py
def __init__(self,
             init_val: Optional[np.ndarray | float] = None,
             timeseries_df: pd.DataFrame = None):
    """
    Args:
        init_val (Optional[np.ndarray | float]):
            starting value(s) at the beginning of the simulation
        timeseries_df (Optional[pd.DataFrame] = None):
            must have columns "date" and "absolute_humidity" --
            "date" entries must correspond to consecutive calendar days
            and must either be strings with `"YYYY-MM-DD"` format or
            `datetime.date` objects -- "value" entries correspond to
            absolute humidity on those days. Identical to
            `FluSubpopSchedules` field of same name.
    """

    super().__init__(init_val)

    self.timeseries_df = timeseries_df

AsympToRecovered

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "IA" to "R" compartment. The functional form is the same across subpopulations.

Source code in CLT_BaseModel/flu_core/flu_components.py
class AsympToRecovered(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "IA" to "R" compartment. The functional form is the same across
    subpopulations.
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.IA_to_R_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.IA_to_R_rate)

BetaReduce

Bases: DynamicVal

"Toy" function representing staged-alert policy that reduces transmission by 50% when more than 5% of the total population is infected. Note: the numbers are completely made up :) The "permanent_lockdown" toggle is to avoid "bang-bang" behavior where the staged-alert policy gets triggered one day and then is off the next, and then is on the day after, and so on... but as the name suggests, it IS permanent. TODO: replace with realistic function.

Source code in CLT_BaseModel/flu_core/flu_components.py
class BetaReduce(clt.DynamicVal):
    """
    "Toy" function representing staged-alert policy
        that reduces transmission by 50% when more than 5%
        of the total population is infected. Note: the
        numbers are completely made up :)
    The "permanent_lockdown" toggle is to avoid "bang-bang"
        behavior where the staged-alert policy gets triggered
        one day and then is off the next, and then is on the
        day after, and so on... but as the name suggests,
        it IS permanent.
    TODO: replace with realistic function.
    """

    def __init__(self, init_val, is_enabled):
        super().__init__(init_val, is_enabled)
        self.permanent_lockdown = False

    def update_current_val(self, state, params):
        if np.sum(sum([state.ISR, state.ISH])) / np.sum(params.total_pop_age_risk) > 0.05:
            self.current_val = .5
            self.permanent_lockdown = True
        else:
            if not self.permanent_lockdown:
                self.current_val = 0.0

DailyVaccines

Bases: Schedule

Source code in CLT_BaseModel/flu_core/flu_components.py
class DailyVaccines(clt.Schedule):

    def __init__(self,
                 init_val: Optional[np.ndarray | float] = None,
                 timeseries_df: pd.DataFrame = None,
                 vax_protection_delay_days: int = 0):
        """
        Args:
            init_val (Optional[np.ndarray | float]):
                starting value(s) at the beginning of the simulation
            timeseries_df (Optional[pd.DataFrame] = None):
                must have "date" and "daily_vaccines" -- "date" entries must
                correspond to consecutive calendar days and must either
                be strings with `"YYYY-MM-DD"` format or `datetime.date`
                objects -- "value" entries correspond to historical
                number vaccinated on those days. Identical to
                `FluSubpopSchedules` field of same name.
            vax_protection_delay_days (int):
                number of days to delay vaccine protection effect.
                Vaccines administered on day X become effective on day X + delay.
        """

        super().__init__(init_val)

        self.timeseries_df = timeseries_df
        self.vax_protection_delay_days = vax_protection_delay_days

    def update_current_val(self, params, current_date: datetime.date) -> None:
        self.current_val = self.timeseries_df.loc[current_date, "daily_vaccines"]

    def postprocess_data_input(self) -> None:
        """
            Converts daily_vaccines column from
            a string representation of a list of lists
            (each day) of format AxR into np.ndarray.
            Shifts dates forward by vax_protection_delay_days
            to model delayed vaccine protection, backfilling
            the beginning with zero entries.
            Pre-indexes the DataFrame by date for O(1) lookups.
        """

        self.timeseries_df['daily_vaccines'] = \
            self.timeseries_df['daily_vaccines'].apply(json.loads)
        self.timeseries_df.loc[:, 'daily_vaccines'] = \
            self.timeseries_df['daily_vaccines'].apply(
                lambda x: np.asarray(x)
                )

        if self.vax_protection_delay_days > 0:
            # Get the original start date and array shape for zero entries
            original_start_date = self.timeseries_df['date'].min()
            zero_array = np.zeros_like(self.timeseries_df['daily_vaccines'].iloc[0])

            # Shift all dates forward by the delay
            self.timeseries_df['date'] = self.timeseries_df['date'].apply(
                lambda d: d + datetime.timedelta(days=self.vax_protection_delay_days)
            )

            # Create backfill rows for the gap at the beginning using pd.date_range
            backfill_dates = pd.date_range(
                start=original_start_date,
                periods=self.vax_protection_delay_days,
                freq='D'
            ).date
            backfill_df = pd.DataFrame({
                'date': backfill_dates,
                'daily_vaccines': [zero_array.copy()] * self.vax_protection_delay_days
            })

            # Concatenate and sort by date
            self.timeseries_df = pd.concat([backfill_df, self.timeseries_df], ignore_index=True)
            self.timeseries_df = self.timeseries_df.sort_values('date').reset_index(drop=True)

        self.timeseries_df = self.timeseries_df.set_index('date')

__init__(init_val: Optional[np.ndarray | float] = None, timeseries_df: pd.DataFrame = None, vax_protection_delay_days: int = 0)

Parameters:

Name Type Description Default
init_val Optional[ndarray | float]

starting value(s) at the beginning of the simulation

None
timeseries_df Optional[pd.DataFrame] = None

must have "date" and "daily_vaccines" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date objects -- "value" entries correspond to historical number vaccinated on those days. Identical to FluSubpopSchedules field of same name.

None
vax_protection_delay_days int

number of days to delay vaccine protection effect. Vaccines administered on day X become effective on day X + delay.

0
Source code in CLT_BaseModel/flu_core/flu_components.py
def __init__(self,
             init_val: Optional[np.ndarray | float] = None,
             timeseries_df: pd.DataFrame = None,
             vax_protection_delay_days: int = 0):
    """
    Args:
        init_val (Optional[np.ndarray | float]):
            starting value(s) at the beginning of the simulation
        timeseries_df (Optional[pd.DataFrame] = None):
            must have "date" and "daily_vaccines" -- "date" entries must
            correspond to consecutive calendar days and must either
            be strings with `"YYYY-MM-DD"` format or `datetime.date`
            objects -- "value" entries correspond to historical
            number vaccinated on those days. Identical to
            `FluSubpopSchedules` field of same name.
        vax_protection_delay_days (int):
            number of days to delay vaccine protection effect.
            Vaccines administered on day X become effective on day X + delay.
    """

    super().__init__(init_val)

    self.timeseries_df = timeseries_df
    self.vax_protection_delay_days = vax_protection_delay_days

postprocess_data_input() -> None

Converts daily_vaccines column from a string representation of a list of lists (each day) of format AxR into np.ndarray. Shifts dates forward by vax_protection_delay_days to model delayed vaccine protection, backfilling the beginning with zero entries. Pre-indexes the DataFrame by date for O(1) lookups.

Source code in CLT_BaseModel/flu_core/flu_components.py
def postprocess_data_input(self) -> None:
    """
        Converts daily_vaccines column from
        a string representation of a list of lists
        (each day) of format AxR into np.ndarray.
        Shifts dates forward by vax_protection_delay_days
        to model delayed vaccine protection, backfilling
        the beginning with zero entries.
        Pre-indexes the DataFrame by date for O(1) lookups.
    """

    self.timeseries_df['daily_vaccines'] = \
        self.timeseries_df['daily_vaccines'].apply(json.loads)
    self.timeseries_df.loc[:, 'daily_vaccines'] = \
        self.timeseries_df['daily_vaccines'].apply(
            lambda x: np.asarray(x)
            )

    if self.vax_protection_delay_days > 0:
        # Get the original start date and array shape for zero entries
        original_start_date = self.timeseries_df['date'].min()
        zero_array = np.zeros_like(self.timeseries_df['daily_vaccines'].iloc[0])

        # Shift all dates forward by the delay
        self.timeseries_df['date'] = self.timeseries_df['date'].apply(
            lambda d: d + datetime.timedelta(days=self.vax_protection_delay_days)
        )

        # Create backfill rows for the gap at the beginning using pd.date_range
        backfill_dates = pd.date_range(
            start=original_start_date,
            periods=self.vax_protection_delay_days,
            freq='D'
        ).date
        backfill_df = pd.DataFrame({
            'date': backfill_dates,
            'daily_vaccines': [zero_array.copy()] * self.vax_protection_delay_days
        })

        # Concatenate and sort by date
        self.timeseries_df = pd.concat([backfill_df, self.timeseries_df], ignore_index=True)
        self.timeseries_df = self.timeseries_df.sort_values('date').reset_index(drop=True)

    self.timeseries_df = self.timeseries_df.set_index('date')

ExposedToAsymp

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "E" to "IA" compartment. The functional form is the same across subpopulations.

Each ExposedToAsymp instance forms a TransitionVariableGroup with a corresponding ExposedToPresymp instance (these two transition variables are jointly distributed).

Source code in CLT_BaseModel/flu_core/flu_components.py
class ExposedToAsymp(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "E" to "IA" compartment. The functional form is the same across
    subpopulations.

    Each ExposedToAsymp instance forms a TransitionVariableGroup with
    a corresponding ExposedToPresymp instance (these two
    transition variables are jointly distributed).
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """
        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.E_to_I_rate * params.E_to_IA_prop)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """
    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.E_to_I_rate * params.E_to_IA_prop)

ExposedToPresymp

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "E" to "IP" compartment. The functional form is the same across subpopulations.

Each ExposedToPresymp instance forms a TransitionVariableGroup with a corresponding ExposedToAsymp instance (these two transition variables are jointly distributed).

Source code in CLT_BaseModel/flu_core/flu_components.py
class ExposedToPresymp(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "E" to "IP" compartment. The functional form is the same across
    subpopulations.

    Each ExposedToPresymp instance forms a TransitionVariableGroup with
    a corresponding ExposedToAsymp instance (these two
    transition variables are jointly distributed).
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.E_to_I_rate * (1 - params.E_to_IA_prop))

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.E_to_I_rate * (1 - params.E_to_IA_prop))

FluContactMatrix

Bases: Schedule

Flu contact matrix.

Attributes:

Name Type Description
timeseries_df DataFrame

must have columns "date", "is_school_day", and "is_work_day" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date object and "is_school_day" and "is_work_day" entries are floats between 0 and 1 indicating if that date is a school day or work day. Identical to FluSubpopSchedules field of same name.

See parent class docstring for other attributes.

Source code in CLT_BaseModel/flu_core/flu_components.py
class FluContactMatrix(clt.Schedule):
    """
    Flu contact matrix.

    Attributes:
        timeseries_df (pd.DataFrame):
            must have columns "date", "is_school_day", and "is_work_day"
            -- "date" entries must correspond to consecutive calendar
            days and must either be strings with `"YYYY-MM-DD"` format
            or `datetime.date` object and "is_school_day" and
            "is_work_day" entries are floats between 0 and 1 indicating if 
            that date is a school day or work day. Identical to 
            `FluSubpopSchedules` field of same name.

    See parent class docstring for other attributes.
    """

    def __init__(self,
                 init_val: Optional[np.ndarray | float] = None,
                 timeseries_df: pd.DataFrame = None):

        super().__init__(init_val)

        self.timeseries_df = timeseries_df

    def update_current_val(self,
                           subpop_params: FluSubpopParams,
                           current_date: datetime.date) -> None:

        try:
            current_row = self.timeseries_df.loc[current_date]
            self.current_val = subpop_params.total_contact_matrix - \
                               (1 - current_row["is_school_day"]) * subpop_params.school_contact_matrix - \
                               (1 - current_row["is_work_day"]) * subpop_params.work_contact_matrix
        except KeyError:
            # print(f"Error: {current_date} is not in `timeseries_df`. Using total contact matrix.")
            self.current_val = subpop_params.total_contact_matrix

    def postprocess_data_input(self) -> None:
        self.timeseries_df = self.timeseries_df.set_index('date')

FluFullMetapopParamsTensors dataclass

Bases: FluTravelParamsTensors

Data container for tensors for FluMetapopModel -- used to store arrays that contain data across all subpopulations (collected from parameters on each location/subpopulation model, as well as from the metapopulation's associated FluMixingParams instance). Note that in contrast to FluTravelParamsTensors, ALL fields in FluSubpopParams are included -- this is for running the simulation via torch.

Attributes:

Name Type Description
num_locations (Tensor, 0 - dimensional)

number of locations (subpopulations) in the metapopulation model and therefore the travel model.

travel_proportions Tensor

L x L array, where L is the number of locations or subpopulations, where element i,j corresponds to the proportion of the population in location i who travels to location j (on average).

See FluSubpopParams docstring for other attributes. Other fields are analogous except they are size (L, A, R) tensors or size 0 tensors.

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluFullMetapopParamsTensors(FluTravelParamsTensors):
    """
    Data container for tensors for `FluMetapopModel` -- used to store arrays that
    contain data across all subpopulations (collected from parameters
    on each location/subpopulation model, as well as from the
    metapopulation's associated `FluMixingParams` instance).
    Note that in contrast to `FluTravelParamsTensors`,
    ALL fields in `FluSubpopParams` are included --
    this is for running the simulation via torch.

    Attributes:
        num_locations (torch.Tensor, 0-dimensional):
            number of locations (subpopulations) in the
            metapopulation model and therefore the travel
            model.
        travel_proportions (torch.Tensor):
            L x L array, where L is the number of locations
            or subpopulations, where element i,j corresponds
            to the proportion of the population in location i
            who travels to location j (on average).

    See `FluSubpopParams` docstring for other attributes.
    Other fields are analogous except they are size (L, A, R)
    tensors or size 0 tensors.
    """

    # non_numerical_params: Optional[dict] = None
    start_real_date: Optional[datetime.date] = None
    beta_baseline: Optional[torch.Tensor] = None
    total_pop_age_risk: Optional[torch.Tensor] = None
    humidity_impact: Optional[torch.Tensor] = None

    inf_induced_saturation: Optional[torch.Tensor] = None
    inf_induced_immune_wane: Optional[torch.Tensor] = None
    vax_induced_saturation: Optional[torch.Tensor] = None
    vax_induced_immune_wane: Optional[torch.Tensor] = None
    inf_induced_inf_risk_reduce: Optional[torch.Tensor] = None
    inf_induced_hosp_risk_reduce: Optional[torch.Tensor] = None
    inf_induced_death_risk_reduce: Optional[torch.Tensor] = None
    vax_induced_inf_risk_reduce: Optional[torch.Tensor] = None
    vax_induced_hosp_risk_reduce: Optional[torch.Tensor] = None
    vax_induced_death_risk_reduce: Optional[torch.Tensor] = None
    vax_induced_inf_risk_reduce_initial: Optional[torch.Tensor] = None
    vax_induced_hosp_risk_reduce_initial: Optional[torch.Tensor] = None
    vax_induced_death_risk_reduce_initial: Optional[torch.Tensor] = None
    adjust_VE_for_seasonal_waning: Optional[torch.Tensor] = True
    vax_protection_delay_days: Optional[torch.Tensor] = 0
    vax_immunity_reset_date_mm_dd: Optional[str] = None

    R_to_S_rate: Optional[torch.Tensor] = None
    E_to_I_rate: Optional[torch.Tensor] = None
    IP_to_IS_rate: Optional[torch.Tensor] = None
    ISR_to_R_rate: Optional[float] = None
    IA_to_R_rate: Optional[float] = None
    ISH_to_H_rate: Optional[float] = None
    HR_to_R_rate: Optional[float] = None
    HD_to_D_rate: Optional[float] = None

    E_to_IA_prop: Optional[torch.Tensor] = None
    IP_to_ISH_prop: Optional[torch.Tensor] = None
    ISH_to_HD_prop: Optional[torch.Tensor] = None

    IP_relative_inf: Optional[torch.Tensor] = None
    IA_relative_inf: Optional[torch.Tensor] = None

    relative_suscept: Optional[torch.Tensor] = None

FluFullMetapopScheduleTensors dataclass

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluFullMetapopScheduleTensors:

    absolute_humidity: Optional[list[torch.tensor]] = None
    is_school_day: Optional[list[torch.tensor]] = None
    is_work_day: Optional[list[torch.tensor]] = None
    daily_vaccines: Optional[list[torch.tensor]] = None
    mobility_modifier: Optional[list[torch.tensor]] = None

FluFullMetapopStateTensors dataclass

Bases: FluTravelStateTensors

Data container for tensors for FluMetapopModel -- used to store arrays that contain data across all subpopulations (collected from each location/subpopulation model). In contrast to FluTravelStateTensors, ALL fields in FluSubpopState are included -- this is for running the simulation via torch.

Attributes:

Name Type Description
flu_contact_matrix torch.Tensor of nonnegative integers

contact matrix for location-age-risk groups -- the lth element holds current_val of FluContactMatrix Schedule for subpopulation l -- this value is a combination of the total contact matrix, the work contact matrix, and the school contact matrix (and the value is adjusted depending on whether the date is a work or school day)

init_vals dict

dictionary of torch.Tensor instances, where keys correspond to "IP", "ISR", "ISH" "IA", "HR", and "HD", and values correspond to their initial values for location-age-risk groups.

See FluSubpopState and FluTravelStateTensors for other attributes -- other attributes here correspond to FluSubpopState, but are size (L, A, R) tensors for location-age-risk or size 0 tensors.

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluFullMetapopStateTensors(FluTravelStateTensors):
    """
    Data container for tensors for `FluMetapopModel` -- used to store arrays that
    contain data across all subpopulations (collected from each
    location/subpopulation model). In contrast to `FluTravelStateTensors`,
    ALL fields in `FluSubpopState` are included -- this is
    for running the simulation via torch.

    Attributes:
        flu_contact_matrix (torch.Tensor of nonnegative integers):
            contact matrix for location-age-risk groups -- the
            lth element holds current_val of `FluContactMatrix`
            `Schedule` for subpopulation l -- this value is a
            combination of the total contact matrix, the
            work contact matrix, and the school contact matrix
            (and the value is adjusted depending on whether
            the date is a work or school day)
        init_vals (dict):
            dictionary of torch.Tensor instances, where keys
            correspond to "IP", "ISR", "ISH" "IA", "HR", and "HD", and values
            correspond to their initial values for location-age-risk
            groups.

    See `FluSubpopState` and `FluTravelStateTensors` for other
        attributes -- other attributes here correspond to
        `FluSubpopState`, but are size (L, A, R) tensors for
        location-age-risk or size 0 tensors.
    """

    # `IP`, `ISR`, `ISH`, `IA`, `HR`, `HD`, `flu_contact_matrix`, `mobility_modifier`
    #   already in parent class
    # Same with `init_vals`

    S: Optional[torch.Tensor] = None
    E: Optional[torch.Tensor] = None
    R: Optional[torch.Tensor] = None
    D: Optional[torch.Tensor] = None

    M: Optional[torch.Tensor] = None
    MV: Optional[torch.Tensor] = None

    absolute_humidity: Optional[float] = None
    daily_vaccines: Optional[torch.Tensor] = None

FluMetapopModel

Bases: MetapopModel, ABC

MetapopModel-derived class specific to flu model.

Source code in CLT_BaseModel/flu_core/flu_components.py
1589
1590
1591
1592
1593
1594
1595
1596
1597
1598
1599
1600
1601
1602
1603
1604
1605
1606
1607
1608
1609
1610
1611
1612
1613
1614
1615
1616
1617
1618
1619
1620
1621
1622
1623
1624
1625
1626
1627
1628
1629
1630
1631
1632
1633
1634
1635
1636
1637
1638
1639
1640
1641
1642
1643
1644
1645
1646
1647
1648
1649
1650
1651
1652
1653
1654
1655
1656
1657
1658
1659
1660
1661
1662
1663
1664
1665
1666
1667
1668
1669
1670
1671
1672
1673
1674
1675
1676
1677
1678
1679
1680
1681
1682
1683
1684
1685
1686
1687
1688
1689
1690
1691
1692
1693
1694
1695
1696
1697
1698
1699
1700
1701
1702
1703
1704
1705
1706
1707
1708
1709
1710
1711
1712
1713
1714
1715
1716
1717
1718
1719
1720
1721
1722
1723
1724
1725
1726
1727
1728
1729
1730
1731
1732
1733
1734
1735
1736
1737
1738
1739
1740
1741
1742
1743
1744
1745
1746
1747
1748
1749
1750
1751
1752
1753
1754
1755
1756
1757
1758
1759
1760
1761
1762
1763
1764
1765
1766
1767
1768
1769
1770
1771
1772
1773
1774
1775
1776
1777
1778
1779
1780
1781
1782
1783
1784
1785
1786
1787
1788
1789
1790
1791
1792
1793
1794
1795
1796
1797
1798
1799
1800
1801
1802
1803
1804
1805
1806
1807
1808
1809
1810
1811
1812
1813
1814
1815
1816
1817
1818
1819
1820
1821
1822
1823
1824
1825
1826
1827
1828
1829
1830
1831
1832
1833
1834
1835
1836
1837
1838
1839
1840
1841
1842
1843
1844
1845
1846
1847
1848
1849
1850
1851
1852
1853
1854
1855
1856
1857
1858
1859
1860
1861
1862
1863
1864
1865
1866
1867
1868
1869
1870
1871
1872
1873
1874
1875
1876
1877
1878
1879
1880
1881
1882
1883
1884
1885
1886
1887
1888
1889
1890
1891
1892
1893
1894
1895
1896
1897
1898
1899
1900
1901
1902
1903
1904
1905
1906
1907
1908
1909
1910
1911
1912
1913
1914
1915
1916
1917
1918
1919
1920
1921
1922
1923
1924
1925
1926
1927
1928
1929
1930
1931
1932
1933
1934
1935
1936
1937
1938
1939
1940
1941
1942
1943
1944
1945
1946
1947
1948
1949
1950
1951
1952
1953
1954
1955
1956
1957
1958
1959
1960
1961
1962
1963
1964
1965
1966
1967
1968
1969
1970
1971
1972
1973
1974
1975
1976
1977
1978
1979
1980
1981
1982
1983
1984
1985
1986
1987
1988
1989
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
2011
2012
2013
2014
2015
2016
2017
2018
2019
2020
2021
2022
2023
2024
2025
2026
2027
2028
2029
2030
2031
2032
2033
2034
2035
2036
2037
2038
2039
2040
2041
2042
2043
2044
2045
2046
2047
2048
2049
2050
2051
2052
2053
2054
2055
2056
2057
2058
2059
2060
2061
2062
2063
2064
2065
2066
2067
2068
2069
2070
2071
2072
2073
2074
2075
2076
2077
2078
2079
2080
2081
2082
class FluMetapopModel(clt.MetapopModel, ABC):
    """
    MetapopModel-derived class specific to flu model.
    """

    def __init__(self,
                 subpop_models: list[dict],
                 mixing_params: FluMixingParams,
                 name: str = ""):

        super().__init__(subpop_models,
                         mixing_params,
                         name)

        # Confirm validity and consistency of `FluMixingParams`
        try:
            num_locations = mixing_params.num_locations
        except KeyError:
            raise FluMetapopModelError("'mixing_params' must contain the key 'num_locations'. \n"
                                       "Please specify it before continuing.")
        if num_locations != len(subpop_models):
            raise FluMetapopModelError("'num_locations' should equal the number of items in \n"
                                       "'subpop_models'. Please amend before continuing.")

        self.travel_state_tensors = FluTravelStateTensors()
        self.update_travel_state_tensors()

        # `FluMixingParams` info is stored on `FluTravelParamsTensors` --
        # this order of operations below is important, because
        # `mixing_params` attribute must be defined before `update_travel_params_tensors()`
        # is called.
        self.mixing_params = mixing_params
        self.travel_params_tensors = FluTravelParamsTensors()
        self.update_travel_params_tensors()

        total_pop_LAR_tensor = self.compute_total_pop_LAR_tensor()

        self.precomputed = FluPrecomputedTensors(total_pop_LAR_tensor,
                                                 self.travel_params_tensors)

        # Generally not used unless using torch version
        self._full_metapop_params_tensors = None
        self._full_metapop_state_tensors = None
        self._full_metapop_schedule_tensors = None

    def check_mobility_input(self) -> None:
        """
        Check that all elements of the mobility matrix are positive,
        and that rows sum to 1.
        """

        travel_proportions = self.mixing_params.travel_proportions
        if np.any(travel_proportions < 0):
            raise FluSubpopModelError(
                f'All entries of the travel matrix must be non-negative:\n {travel_proportions}.')

        if not(np.allclose(travel_proportions.sum(axis=1), 1)):
            raise FluSubpopModelError(
                f'Rows of the travel matrix must all sum to 1:\n {travel_proportions}.')

    def run_input_checks(self) -> None:
        """
        Check the following:
            - rows of mobility matrix must sum to 1
            - mobility matrix entries are between 0 and 1
        """

        # Mobility matrix check
        self.check_mobility_input()

    def modify_subpop_params(self,
                             subpop_name: str,
                             updates_dict: dict):
        """
        This method lets users safely modify a single subpopulation
        parameters field; the metapopulation-wide tensors are updated
        automatically afterward.

        In a `FluMetapopModel`, subpopulation parameters are combined into
        (L, A, R) tensors across L subpopulations.`FluSubpopParams` is a frozen
        dataclass to avoid users naively changing parameter values and getting
        undesirable results -- thus, `FluSubpopParams` on a subpopulation
        model cannot be updated directly.

        Parameters:
            subpop_name (str):
               Value must match the `name` attribute of one of the
               `FluSubpopModel` instances contained in this metapopulation
                model's `subpop_models` attribute.
            updates_dict (dict):
                Dictionary specifying values to update in a
                `FluSubpopParams` instance -- keys must match the
                field names of `FluSubpopParams`.
        """

        # Since `FluSubpopParams` is frozen, we return a new instance
        #   with the reflected updates
        self.subpop_models[subpop_name].params = clt.updated_dataclass(
            self.subpop_models[subpop_name].params, updates_dict
        )

        self.update_travel_params_tensors()

        # Adding this for extra safety in case the user does not
        # call `get_flu_torch_inputs` for accessing the
        # `FullMetapopParams` instance.

        # If this attribute is not `None`, it means we are using
        # the `torch` implementation, and we should update the
        # corresponding `FullMetapopParams` instance with the new
        # `FluMixingParams` values.
        if self._full_metapop_params_tensors:
            self.update_full_metapop_params_tensors()

    def modify_mixing_params(self,
                             updates_dict: dict):
        """
        This method lets users safely modify flu mixing parameters;
        the metapopulation-wide tensors are updated automatically afterward.
        `FluMixingParams` is a frozen dataclass to avoid users
        naively changing parameter values and getting undesirable results --
        thus, `FluMixingParams` cannot be updated directly.

        Parameters:
            updates_dict (dict):
                Dictionary specifying values to update in a
                `FluSubpopParams` instance -- keys must match the
                field names of `FluSubpopParams`. 
        """

        self.mixing_params = clt.updated_dataclass(self.mixing_params, updates_dict)
        self.update_travel_params_tensors()

        nonlocal_travel_prop = self.travel_params_tensors.travel_proportions.clone().fill_diagonal_(0.0)

        self.precomputed.sum_residents_nonlocal_travel_prop = nonlocal_travel_prop.sum(dim=1)

        # Adding this for extra safety in case the user does not
        # call `get_flu_torch_inputs` for accessing the
        # `FullMetapopParams` instance.

        # If this attribute is not `None`, it means we are using
        # the `torch` implementation, and we should update the
        # corresponding `FullMetapopParams` instance with the new
        # `FluMixingParams` values.
        if self._full_metapop_params_tensors:
            self.update_full_metapop_params_tensors()

    def compute_total_pop_LAR_tensor(self) -> torch.tensor:
        """
        For each subpopulation, sum initial values of population
        in each compartment for age-risk groups. Store all information
        as tensor and return tensor.

        Returns:
        --------
        torch.tensor of size (L, A, R):
            Total population (across all compartments) for
            location-age-risk (l, a, r).
        """

        # ORDER MATTERS! USE ORDERED DICTIONARY HERE
        #   to preserve correct index order in tensors!
        #   See `update_travel_params_tensors` for detailed note.
        subpop_models_ordered = self._subpop_models_ordered

        total_pop_LAR_tensor = torch.zeros(self.travel_params_tensors.num_locations,
                                           self.travel_params_tensors.num_age_groups,
                                           self.travel_params_tensors.num_risk_groups)

        # All subpop models should have the same compartments' keys
        for name in subpop_models_ordered[0].compartments.keys():

            metapop_vals = []

            for model in subpop_models_ordered.values():
                compartment = getattr(model.compartments, name)
                metapop_vals.append(compartment.current_val)

            total_pop_LAR_tensor = total_pop_LAR_tensor + torch.tensor(np.asarray(metapop_vals))

        return total_pop_LAR_tensor

    def update_state_tensors(self,
                             target: FluTravelStateTensors) -> None:
        """
        Update `target` instance in-place with current simulation
        values. Each field of `target` corresponds to a field in
        `FluSubpopState`, and contains either a tensor of size
        (L, A, R) or a tensor of size (L), where (l, a, r) refers to
        location-age-risk.
        """

        # ORDER MATTERS! USE ORDERED DICTIONARY HERE
        #   to preserve correct index order in tensors!
        #   See `update_travel_params_tensors` for detailed note.
        subpop_models_ordered = self._subpop_models_ordered

        for field in fields(target):

            name = field.name

            # FluTravelStateTensors has an attribute
            #   that is a dictionary called `init_vals` --
            #   disregard, as this only used to store
            #   initial values for resetting, but is not
            #   used in the travel model computation
            if name == "init_vals":
                continue

            metapop_vals = []

            for model in subpop_models_ordered.values():
                current_val = getattr(model.state, name)
                metapop_vals.append(current_val)

            # Probably want to update this to be cleaner...
            # `SubpopState` fields that correspond to `Schedule` instances
            # have initial values of `None` -- but we cannot build a tensor
            # with `None` values, so we convert values to 0s.
            if any(v is None for v in metapop_vals):
                setattr(target, name, torch.tensor(np.full(np.shape(metapop_vals), 0.0)))
            else:
                setattr(target, name, torch.tensor(np.asarray(metapop_vals)))

            # Only fields corresponding to `Schedule` instances can be
            # size (L) -- this is because the schedule value may be scalar for
            # each subpopulation. Other fields should all be size (L, A, R). 

    def update_travel_state_tensors(self) -> None:
        """
        Update `travel_state_tensors` attribute in-place.
        `FluTravelStateTensors` only has fields corresponding
        to state variables relevant for the travel model.
        Converts subpopulation-specific state to
        tensors of size (L, A, R) for location-age-risk
        (except for a few exceptions that have different dimensions).
        """

        self.update_state_tensors(self.travel_state_tensors)

    def update_full_metapop_state_tensors(self) -> None:
        """
        Update `_full_metapop_state_tensors` attribute in-place.
        `FluFullMetapopStateTensors` has fields corresponding
        to all state variables in the simulation.
        Converts subpopulation-specific state to
        tensors of size (L, A, R) for location-age-risk
        (except for a few exceptions that have different dimensions).
        """

        if self._full_metapop_state_tensors is None:
            self._full_metapop_state_tensors = FluFullMetapopStateTensors()
        self.update_state_tensors(self._full_metapop_state_tensors)

    def update_params_tensors(self,
                              target: FluTravelParamsTensors) -> FluTravelParamsTensors:
        """
        Update `target` in-place. Converts subpopulation-specific
        parameters to tensors of size (L, A, R) for location-age-risk,
        except for `num_locations` and `travel_proportions`, which
        have size 1 and (L, L) respectively.
        """

        # USE THE ORDERED DICTIONARY HERE FOR SAFETY!
        #   AGAIN, ORDER MATTERS BECAUSE ORDER DETERMINES
        #   THE SUBPOPULATION INDEX IN THE METAPOPULATION
        #   TENSOR!
        subpop_models_ordered = self._subpop_models_ordered

        # Subpop models should have the same A, R so grab
        #   from the first subpop model
        A = subpop_models_ordered[0].params.num_age_groups
        R = subpop_models_ordered[0].params.num_risk_groups

        for field in fields(target):

            name = field.name
            is_non_numerical = False

            metapop_vals = []

            if name == "num_locations" or name == "travel_proportions":
                setattr(target, name, torch.tensor(getattr(self.mixing_params, name)))

            else:

                for model in subpop_models_ordered.values():
                    metapop_vals.append(getattr(model.params, name))

                # If all values are equal to each other, then
                #   simply store the first value (since its value is common
                #   across metapopulations)
                first_val = metapop_vals[0]
                if isinstance(first_val, str) or isinstance(first_val, datetime.date):
                    is_non_numerical = True
                    if all(x == first_val for x in metapop_vals):
                        metapop_vals = first_val
                    else:
                        raise FluMetapopModelError(
                            f"Error: non-numerical parameter '{name}' has values that differ "
                            "across subpopulations; values should be the same."
                        )
                else:
                    if all(np.allclose(x, first_val) for x in metapop_vals):
                        metapop_vals = first_val

                # Converting list of arrays to tensors is slow --
                #   better to convert to array first
                if isinstance(metapop_vals, list):
                    metapop_vals = np.asarray(metapop_vals)
                    # metapop_vals = np.stack([clt.to_AR_array(x, A, R) for x in metapop_vals])

                if is_non_numerical:
                    setattr(target, name, metapop_vals)
                else:
                    setattr(target, name, torch.tensor(metapop_vals))

        # Convert all tensors to correct size!
        target.standardize_shapes()

    def update_travel_params_tensors(self) -> None:
        """
        Update `travel_params_tensors` attribute in-place.
        `FluTravelParamsTensors` only has fields corresponding
        to parameters relevant for the travel model.
        Converts subpopulation-specific parameters to
        tensors of size (L, A, R) for location-age-risk
        (except for a few exceptions that have different dimensions).
        """

        self.update_params_tensors(target=self.travel_params_tensors)

    def update_full_metapop_params_tensors(self) -> None:
        """
        Update `_full_metapop_params_tensors` attribute in-place.
        `FluFullMetapopParamsTensors` has fields corresponding
        to all parameters in the simulation. Converts subpopulation-specific
        parameters to tensors of size (L, A, R) for location-age-risk
        (except for a few exceptions that have different dimensions).
        """

        if self._full_metapop_params_tensors is None:
            self._full_metapop_params_tensors = FluFullMetapopParamsTensors()
        self.update_params_tensors(target=self._full_metapop_params_tensors)

    def apply_inter_subpop_updates(self) -> None:
        """
        Update the `FluTravelStateTensors` according to the simulation state
        and compute the total mixing exposure, which includes across-subpopulation
        mixing/travel. Update the `total_mixing_exposure` attribute on each
        subpopulation's `SusceptibleToExposed` instance accordingly, so each
        of these transition variables can compute its transition rate.

        See `apply_inter_subpop_updates` on `MetapopModel` base class
        for logic of how/when this is called in the simulation.
        """

        self.update_travel_state_tensors()

        total_mixing_exposure = compute_total_mixing_exposure(self.travel_state_tensors,
                                                              self.travel_params_tensors,
                                                              self.precomputed)

        # Again, `self.subpop_models` is an ordered dictionary --
        #   so iterating over the dictionary like this is well-defined
        #   and responsible -- the order is important because it
        #   determines the order (index) in any metapopulation tensors
        subpop_models = self._subpop_models_ordered

        # Updates `total_mixing_exposure` attribute on each `SusceptibleToExposed`
        # instance -- this value captures across-population travel/mixing.
        for i in range(len(subpop_models)):
            subpop_models.values()[i].transition_variables.S_to_E.total_mixing_exposure = \
                total_mixing_exposure[i, :, :]

    def setup_full_metapop_schedule_tensors(self):
        """
        Creates `FluFullMetapopScheduleTensors` instance and assigns to
        `_full_metapop_schedule_tensors` attribute.

        For the metapopulation model's L locations/subpopulations, for each day,
        each value-related column in each schedule is either a float or
        array of size (A, R) for age-risk groups.

        We aggregate and reformat this schedule information and put it
        into a `FluFullMetapopScheduleTensors` instance, where fields
        correspond to a schedule value, and values are lists of tensors of
        size (L, A, R). The ith element of each list corresponds to the
        ith simulation day.
        """

        self._full_metapop_schedule_tensors = FluFullMetapopScheduleTensors()

        L = self.precomputed.L
        A = self.precomputed.A
        R = self.precomputed.R

        # Note: there is probably a more consistent way to do this,
        # because now `flu_contact_matrix` has two values: "is_school_day"
        # and "is_work_day" -- other schedules' dataframes only have one
        # relevant column value rather than two
        for item in [("absolute_humidity", "absolute_humidity"),
                     ("flu_contact_matrix", "is_school_day"),
                     ("flu_contact_matrix", "is_work_day"),
                     ("daily_vaccines", "daily_vaccines"),
                     ("mobility_modifier", "mobility_modifier")]:

            schedule_name = item[0]
            values_column_name = item[1]

            metapop_vals = []

            for subpop_model in self._subpop_models_ordered.values():
                df = subpop_model.schedules[schedule_name].timeseries_df

                # Using the `start_real_date` specification given in subpop's `SimulationSettings`,
                # extract the relevant part of the dataframe with dates >= the simulation start date.
                # Note that `start_real_date` should be the same for each subpopulation
                start_date = datetime.datetime.strptime(subpop_model.simulation_settings.start_real_date, "%Y-%m-%d")

                # If schedule uses day_of_week scheduling, we need to create the full date range
                # for the schedule dataframe
                if subpop_model.schedules[schedule_name].is_day_of_week_schedule:
                    df = create_timeseries_df_from_day_of_week_schedule(
                        df, start_date)

                df["simulation_day"] = (pd.to_datetime(df.index, format="%Y-%m-%d") - start_date).to_series().dt.days.values
                df = df[df["simulation_day"] >= 0]

                # Make each day's value an A x R array
                # Pandas complains about `SettingWithCopyWarning` so we work on a copy explicitly to stop it
                #   from complaining...
                df = df.copy()

                if schedule_name in ['daily_vaccines', 'mobility_modifier']:
                    # daily_vaccines and mobility_modifier are already given as A x R arrays
                    if df[values_column_name].values[0].shape != (A, R):
                        raise ValueError(f"Error: {schedule_name} arrays must have shape ({A}, {R}). " \
                            f"Current input has shape {df[values_column_name].values[0].shape}.")
                else:
                    df[values_column_name] = df[values_column_name].astype(object)
                    df.loc[:, values_column_name] = df[values_column_name].apply(
                        lambda x, A=A, R=R: np.broadcast_to(np.asarray(x).reshape(1, 1), (A, R))
                    )

                metapop_vals.append(np.asarray(df[values_column_name]))

            # IMPORTANT: tedious array/tensor shape/size manipulation here
            # metapop_vals: list of L arrays, each shape (num_days, A, R)
            # We need to transpose this... to be a list of num_days tensors, of size L x A x R
            num_items = metapop_vals[0].shape[0]

            # This is ugly and inefficient -- but at least we only do this once, when we get the initial
            #   state of a metapopulation model in tensor form
            transposed_metapop_vals = [torch.tensor(np.array([metapop_vals[l][i] for l in range(L)])) for i in
                                       range(num_items)]

            setattr(self._full_metapop_schedule_tensors, values_column_name, transposed_metapop_vals)

    def get_flu_torch_inputs(self) -> dict:
        """
        Prepares and returns metapopulation simulation data in tensor format
        that can be directly used for `torch` implementation.

        Returns:
             d (dict):
                Has keys "state_tensors", "params_tensors", "schedule_tensors",
                and "precomputed". Corresponds to `FluFullMetapopStateTensors`,
                `FluFullMetapopParamsTensors`, `FluFullMetapopScheduleTensors`,
                and `FluPrecomputedTensors` instances respectively.
        """

        # Note: does not support dynamic variables (yet). If want to
        #   run pytorch with dynamic variables, will need to create
        #   a method similar to `setup_full_metapop_schedule_tensors`
        #   but for dynamic variables. Also note that we cannot differentiate
        #   with respect to dynamic variables that are discontinuous
        #   (e.g. a 0-1 intervention) -- so we cannot optimize discontinuous
        #   dynamic variables.

        self.update_full_metapop_state_tensors()
        self.update_full_metapop_params_tensors()
        self._full_metapop_params_tensors.standardize_shapes()
        self.setup_full_metapop_schedule_tensors()

        d = {}

        d["state_tensors"] = copy.deepcopy(self._full_metapop_state_tensors)
        d["params_tensors"] = copy.deepcopy(self._full_metapop_params_tensors)
        d["schedule_tensors"] = copy.deepcopy(self._full_metapop_schedule_tensors)
        d["precomputed"] = copy.deepcopy(self.precomputed)

        return d

apply_inter_subpop_updates() -> None

Update the FluTravelStateTensors according to the simulation state and compute the total mixing exposure, which includes across-subpopulation mixing/travel. Update the total_mixing_exposure attribute on each subpopulation's SusceptibleToExposed instance accordingly, so each of these transition variables can compute its transition rate.

See apply_inter_subpop_updates on MetapopModel base class for logic of how/when this is called in the simulation.

Source code in CLT_BaseModel/flu_core/flu_components.py
def apply_inter_subpop_updates(self) -> None:
    """
    Update the `FluTravelStateTensors` according to the simulation state
    and compute the total mixing exposure, which includes across-subpopulation
    mixing/travel. Update the `total_mixing_exposure` attribute on each
    subpopulation's `SusceptibleToExposed` instance accordingly, so each
    of these transition variables can compute its transition rate.

    See `apply_inter_subpop_updates` on `MetapopModel` base class
    for logic of how/when this is called in the simulation.
    """

    self.update_travel_state_tensors()

    total_mixing_exposure = compute_total_mixing_exposure(self.travel_state_tensors,
                                                          self.travel_params_tensors,
                                                          self.precomputed)

    # Again, `self.subpop_models` is an ordered dictionary --
    #   so iterating over the dictionary like this is well-defined
    #   and responsible -- the order is important because it
    #   determines the order (index) in any metapopulation tensors
    subpop_models = self._subpop_models_ordered

    # Updates `total_mixing_exposure` attribute on each `SusceptibleToExposed`
    # instance -- this value captures across-population travel/mixing.
    for i in range(len(subpop_models)):
        subpop_models.values()[i].transition_variables.S_to_E.total_mixing_exposure = \
            total_mixing_exposure[i, :, :]

check_mobility_input() -> None

Check that all elements of the mobility matrix are positive, and that rows sum to 1.

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_mobility_input(self) -> None:
    """
    Check that all elements of the mobility matrix are positive,
    and that rows sum to 1.
    """

    travel_proportions = self.mixing_params.travel_proportions
    if np.any(travel_proportions < 0):
        raise FluSubpopModelError(
            f'All entries of the travel matrix must be non-negative:\n {travel_proportions}.')

    if not(np.allclose(travel_proportions.sum(axis=1), 1)):
        raise FluSubpopModelError(
            f'Rows of the travel matrix must all sum to 1:\n {travel_proportions}.')

compute_total_pop_LAR_tensor() -> torch.tensor

For each subpopulation, sum initial values of population in each compartment for age-risk groups. Store all information as tensor and return tensor.

Returns:

torch.tensor of size (L, A, R): Total population (across all compartments) for location-age-risk (l, a, r).

Source code in CLT_BaseModel/flu_core/flu_components.py
def compute_total_pop_LAR_tensor(self) -> torch.tensor:
    """
    For each subpopulation, sum initial values of population
    in each compartment for age-risk groups. Store all information
    as tensor and return tensor.

    Returns:
    --------
    torch.tensor of size (L, A, R):
        Total population (across all compartments) for
        location-age-risk (l, a, r).
    """

    # ORDER MATTERS! USE ORDERED DICTIONARY HERE
    #   to preserve correct index order in tensors!
    #   See `update_travel_params_tensors` for detailed note.
    subpop_models_ordered = self._subpop_models_ordered

    total_pop_LAR_tensor = torch.zeros(self.travel_params_tensors.num_locations,
                                       self.travel_params_tensors.num_age_groups,
                                       self.travel_params_tensors.num_risk_groups)

    # All subpop models should have the same compartments' keys
    for name in subpop_models_ordered[0].compartments.keys():

        metapop_vals = []

        for model in subpop_models_ordered.values():
            compartment = getattr(model.compartments, name)
            metapop_vals.append(compartment.current_val)

        total_pop_LAR_tensor = total_pop_LAR_tensor + torch.tensor(np.asarray(metapop_vals))

    return total_pop_LAR_tensor

get_flu_torch_inputs() -> dict

Prepares and returns metapopulation simulation data in tensor format that can be directly used for torch implementation.

Returns:

Name Type Description
d dict

Has keys "state_tensors", "params_tensors", "schedule_tensors", and "precomputed". Corresponds to FluFullMetapopStateTensors, FluFullMetapopParamsTensors, FluFullMetapopScheduleTensors, and FluPrecomputedTensors instances respectively.

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_flu_torch_inputs(self) -> dict:
    """
    Prepares and returns metapopulation simulation data in tensor format
    that can be directly used for `torch` implementation.

    Returns:
         d (dict):
            Has keys "state_tensors", "params_tensors", "schedule_tensors",
            and "precomputed". Corresponds to `FluFullMetapopStateTensors`,
            `FluFullMetapopParamsTensors`, `FluFullMetapopScheduleTensors`,
            and `FluPrecomputedTensors` instances respectively.
    """

    # Note: does not support dynamic variables (yet). If want to
    #   run pytorch with dynamic variables, will need to create
    #   a method similar to `setup_full_metapop_schedule_tensors`
    #   but for dynamic variables. Also note that we cannot differentiate
    #   with respect to dynamic variables that are discontinuous
    #   (e.g. a 0-1 intervention) -- so we cannot optimize discontinuous
    #   dynamic variables.

    self.update_full_metapop_state_tensors()
    self.update_full_metapop_params_tensors()
    self._full_metapop_params_tensors.standardize_shapes()
    self.setup_full_metapop_schedule_tensors()

    d = {}

    d["state_tensors"] = copy.deepcopy(self._full_metapop_state_tensors)
    d["params_tensors"] = copy.deepcopy(self._full_metapop_params_tensors)
    d["schedule_tensors"] = copy.deepcopy(self._full_metapop_schedule_tensors)
    d["precomputed"] = copy.deepcopy(self.precomputed)

    return d

modify_mixing_params(updates_dict: dict)

This method lets users safely modify flu mixing parameters; the metapopulation-wide tensors are updated automatically afterward. FluMixingParams is a frozen dataclass to avoid users naively changing parameter values and getting undesirable results -- thus, FluMixingParams cannot be updated directly.

Parameters:

Name Type Description Default
updates_dict dict

Dictionary specifying values to update in a FluSubpopParams instance -- keys must match the field names of FluSubpopParams.

required
Source code in CLT_BaseModel/flu_core/flu_components.py
def modify_mixing_params(self,
                         updates_dict: dict):
    """
    This method lets users safely modify flu mixing parameters;
    the metapopulation-wide tensors are updated automatically afterward.
    `FluMixingParams` is a frozen dataclass to avoid users
    naively changing parameter values and getting undesirable results --
    thus, `FluMixingParams` cannot be updated directly.

    Parameters:
        updates_dict (dict):
            Dictionary specifying values to update in a
            `FluSubpopParams` instance -- keys must match the
            field names of `FluSubpopParams`. 
    """

    self.mixing_params = clt.updated_dataclass(self.mixing_params, updates_dict)
    self.update_travel_params_tensors()

    nonlocal_travel_prop = self.travel_params_tensors.travel_proportions.clone().fill_diagonal_(0.0)

    self.precomputed.sum_residents_nonlocal_travel_prop = nonlocal_travel_prop.sum(dim=1)

    # Adding this for extra safety in case the user does not
    # call `get_flu_torch_inputs` for accessing the
    # `FullMetapopParams` instance.

    # If this attribute is not `None`, it means we are using
    # the `torch` implementation, and we should update the
    # corresponding `FullMetapopParams` instance with the new
    # `FluMixingParams` values.
    if self._full_metapop_params_tensors:
        self.update_full_metapop_params_tensors()

modify_subpop_params(subpop_name: str, updates_dict: dict)

This method lets users safely modify a single subpopulation parameters field; the metapopulation-wide tensors are updated automatically afterward.

In a FluMetapopModel, subpopulation parameters are combined into (L, A, R) tensors across L subpopulations.FluSubpopParams is a frozen dataclass to avoid users naively changing parameter values and getting undesirable results -- thus, FluSubpopParams on a subpopulation model cannot be updated directly.

Parameters:

Name Type Description Default
subpop_name str

Value must match the name attribute of one of the FluSubpopModel instances contained in this metapopulation model's subpop_models attribute.

required
updates_dict dict

Dictionary specifying values to update in a FluSubpopParams instance -- keys must match the field names of FluSubpopParams.

required
Source code in CLT_BaseModel/flu_core/flu_components.py
def modify_subpop_params(self,
                         subpop_name: str,
                         updates_dict: dict):
    """
    This method lets users safely modify a single subpopulation
    parameters field; the metapopulation-wide tensors are updated
    automatically afterward.

    In a `FluMetapopModel`, subpopulation parameters are combined into
    (L, A, R) tensors across L subpopulations.`FluSubpopParams` is a frozen
    dataclass to avoid users naively changing parameter values and getting
    undesirable results -- thus, `FluSubpopParams` on a subpopulation
    model cannot be updated directly.

    Parameters:
        subpop_name (str):
           Value must match the `name` attribute of one of the
           `FluSubpopModel` instances contained in this metapopulation
            model's `subpop_models` attribute.
        updates_dict (dict):
            Dictionary specifying values to update in a
            `FluSubpopParams` instance -- keys must match the
            field names of `FluSubpopParams`.
    """

    # Since `FluSubpopParams` is frozen, we return a new instance
    #   with the reflected updates
    self.subpop_models[subpop_name].params = clt.updated_dataclass(
        self.subpop_models[subpop_name].params, updates_dict
    )

    self.update_travel_params_tensors()

    # Adding this for extra safety in case the user does not
    # call `get_flu_torch_inputs` for accessing the
    # `FullMetapopParams` instance.

    # If this attribute is not `None`, it means we are using
    # the `torch` implementation, and we should update the
    # corresponding `FullMetapopParams` instance with the new
    # `FluMixingParams` values.
    if self._full_metapop_params_tensors:
        self.update_full_metapop_params_tensors()

run_input_checks() -> None

Check the following
  • rows of mobility matrix must sum to 1
  • mobility matrix entries are between 0 and 1
Source code in CLT_BaseModel/flu_core/flu_components.py
def run_input_checks(self) -> None:
    """
    Check the following:
        - rows of mobility matrix must sum to 1
        - mobility matrix entries are between 0 and 1
    """

    # Mobility matrix check
    self.check_mobility_input()

setup_full_metapop_schedule_tensors()

Creates FluFullMetapopScheduleTensors instance and assigns to _full_metapop_schedule_tensors attribute.

For the metapopulation model's L locations/subpopulations, for each day, each value-related column in each schedule is either a float or array of size (A, R) for age-risk groups.

We aggregate and reformat this schedule information and put it into a FluFullMetapopScheduleTensors instance, where fields correspond to a schedule value, and values are lists of tensors of size (L, A, R). The ith element of each list corresponds to the ith simulation day.

Source code in CLT_BaseModel/flu_core/flu_components.py
def setup_full_metapop_schedule_tensors(self):
    """
    Creates `FluFullMetapopScheduleTensors` instance and assigns to
    `_full_metapop_schedule_tensors` attribute.

    For the metapopulation model's L locations/subpopulations, for each day,
    each value-related column in each schedule is either a float or
    array of size (A, R) for age-risk groups.

    We aggregate and reformat this schedule information and put it
    into a `FluFullMetapopScheduleTensors` instance, where fields
    correspond to a schedule value, and values are lists of tensors of
    size (L, A, R). The ith element of each list corresponds to the
    ith simulation day.
    """

    self._full_metapop_schedule_tensors = FluFullMetapopScheduleTensors()

    L = self.precomputed.L
    A = self.precomputed.A
    R = self.precomputed.R

    # Note: there is probably a more consistent way to do this,
    # because now `flu_contact_matrix` has two values: "is_school_day"
    # and "is_work_day" -- other schedules' dataframes only have one
    # relevant column value rather than two
    for item in [("absolute_humidity", "absolute_humidity"),
                 ("flu_contact_matrix", "is_school_day"),
                 ("flu_contact_matrix", "is_work_day"),
                 ("daily_vaccines", "daily_vaccines"),
                 ("mobility_modifier", "mobility_modifier")]:

        schedule_name = item[0]
        values_column_name = item[1]

        metapop_vals = []

        for subpop_model in self._subpop_models_ordered.values():
            df = subpop_model.schedules[schedule_name].timeseries_df

            # Using the `start_real_date` specification given in subpop's `SimulationSettings`,
            # extract the relevant part of the dataframe with dates >= the simulation start date.
            # Note that `start_real_date` should be the same for each subpopulation
            start_date = datetime.datetime.strptime(subpop_model.simulation_settings.start_real_date, "%Y-%m-%d")

            # If schedule uses day_of_week scheduling, we need to create the full date range
            # for the schedule dataframe
            if subpop_model.schedules[schedule_name].is_day_of_week_schedule:
                df = create_timeseries_df_from_day_of_week_schedule(
                    df, start_date)

            df["simulation_day"] = (pd.to_datetime(df.index, format="%Y-%m-%d") - start_date).to_series().dt.days.values
            df = df[df["simulation_day"] >= 0]

            # Make each day's value an A x R array
            # Pandas complains about `SettingWithCopyWarning` so we work on a copy explicitly to stop it
            #   from complaining...
            df = df.copy()

            if schedule_name in ['daily_vaccines', 'mobility_modifier']:
                # daily_vaccines and mobility_modifier are already given as A x R arrays
                if df[values_column_name].values[0].shape != (A, R):
                    raise ValueError(f"Error: {schedule_name} arrays must have shape ({A}, {R}). " \
                        f"Current input has shape {df[values_column_name].values[0].shape}.")
            else:
                df[values_column_name] = df[values_column_name].astype(object)
                df.loc[:, values_column_name] = df[values_column_name].apply(
                    lambda x, A=A, R=R: np.broadcast_to(np.asarray(x).reshape(1, 1), (A, R))
                )

            metapop_vals.append(np.asarray(df[values_column_name]))

        # IMPORTANT: tedious array/tensor shape/size manipulation here
        # metapop_vals: list of L arrays, each shape (num_days, A, R)
        # We need to transpose this... to be a list of num_days tensors, of size L x A x R
        num_items = metapop_vals[0].shape[0]

        # This is ugly and inefficient -- but at least we only do this once, when we get the initial
        #   state of a metapopulation model in tensor form
        transposed_metapop_vals = [torch.tensor(np.array([metapop_vals[l][i] for l in range(L)])) for i in
                                   range(num_items)]

        setattr(self._full_metapop_schedule_tensors, values_column_name, transposed_metapop_vals)

update_full_metapop_params_tensors() -> None

Update _full_metapop_params_tensors attribute in-place. FluFullMetapopParamsTensors has fields corresponding to all parameters in the simulation. Converts subpopulation-specific parameters to tensors of size (L, A, R) for location-age-risk (except for a few exceptions that have different dimensions).

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_full_metapop_params_tensors(self) -> None:
    """
    Update `_full_metapop_params_tensors` attribute in-place.
    `FluFullMetapopParamsTensors` has fields corresponding
    to all parameters in the simulation. Converts subpopulation-specific
    parameters to tensors of size (L, A, R) for location-age-risk
    (except for a few exceptions that have different dimensions).
    """

    if self._full_metapop_params_tensors is None:
        self._full_metapop_params_tensors = FluFullMetapopParamsTensors()
    self.update_params_tensors(target=self._full_metapop_params_tensors)

update_full_metapop_state_tensors() -> None

Update _full_metapop_state_tensors attribute in-place. FluFullMetapopStateTensors has fields corresponding to all state variables in the simulation. Converts subpopulation-specific state to tensors of size (L, A, R) for location-age-risk (except for a few exceptions that have different dimensions).

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_full_metapop_state_tensors(self) -> None:
    """
    Update `_full_metapop_state_tensors` attribute in-place.
    `FluFullMetapopStateTensors` has fields corresponding
    to all state variables in the simulation.
    Converts subpopulation-specific state to
    tensors of size (L, A, R) for location-age-risk
    (except for a few exceptions that have different dimensions).
    """

    if self._full_metapop_state_tensors is None:
        self._full_metapop_state_tensors = FluFullMetapopStateTensors()
    self.update_state_tensors(self._full_metapop_state_tensors)

update_params_tensors(target: FluTravelParamsTensors) -> FluTravelParamsTensors

Update target in-place. Converts subpopulation-specific parameters to tensors of size (L, A, R) for location-age-risk, except for num_locations and travel_proportions, which have size 1 and (L, L) respectively.

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_params_tensors(self,
                          target: FluTravelParamsTensors) -> FluTravelParamsTensors:
    """
    Update `target` in-place. Converts subpopulation-specific
    parameters to tensors of size (L, A, R) for location-age-risk,
    except for `num_locations` and `travel_proportions`, which
    have size 1 and (L, L) respectively.
    """

    # USE THE ORDERED DICTIONARY HERE FOR SAFETY!
    #   AGAIN, ORDER MATTERS BECAUSE ORDER DETERMINES
    #   THE SUBPOPULATION INDEX IN THE METAPOPULATION
    #   TENSOR!
    subpop_models_ordered = self._subpop_models_ordered

    # Subpop models should have the same A, R so grab
    #   from the first subpop model
    A = subpop_models_ordered[0].params.num_age_groups
    R = subpop_models_ordered[0].params.num_risk_groups

    for field in fields(target):

        name = field.name
        is_non_numerical = False

        metapop_vals = []

        if name == "num_locations" or name == "travel_proportions":
            setattr(target, name, torch.tensor(getattr(self.mixing_params, name)))

        else:

            for model in subpop_models_ordered.values():
                metapop_vals.append(getattr(model.params, name))

            # If all values are equal to each other, then
            #   simply store the first value (since its value is common
            #   across metapopulations)
            first_val = metapop_vals[0]
            if isinstance(first_val, str) or isinstance(first_val, datetime.date):
                is_non_numerical = True
                if all(x == first_val for x in metapop_vals):
                    metapop_vals = first_val
                else:
                    raise FluMetapopModelError(
                        f"Error: non-numerical parameter '{name}' has values that differ "
                        "across subpopulations; values should be the same."
                    )
            else:
                if all(np.allclose(x, first_val) for x in metapop_vals):
                    metapop_vals = first_val

            # Converting list of arrays to tensors is slow --
            #   better to convert to array first
            if isinstance(metapop_vals, list):
                metapop_vals = np.asarray(metapop_vals)
                # metapop_vals = np.stack([clt.to_AR_array(x, A, R) for x in metapop_vals])

            if is_non_numerical:
                setattr(target, name, metapop_vals)
            else:
                setattr(target, name, torch.tensor(metapop_vals))

    # Convert all tensors to correct size!
    target.standardize_shapes()

update_state_tensors(target: FluTravelStateTensors) -> None

Update target instance in-place with current simulation values. Each field of target corresponds to a field in FluSubpopState, and contains either a tensor of size (L, A, R) or a tensor of size (L), where (l, a, r) refers to location-age-risk.

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_state_tensors(self,
                         target: FluTravelStateTensors) -> None:
    """
    Update `target` instance in-place with current simulation
    values. Each field of `target` corresponds to a field in
    `FluSubpopState`, and contains either a tensor of size
    (L, A, R) or a tensor of size (L), where (l, a, r) refers to
    location-age-risk.
    """

    # ORDER MATTERS! USE ORDERED DICTIONARY HERE
    #   to preserve correct index order in tensors!
    #   See `update_travel_params_tensors` for detailed note.
    subpop_models_ordered = self._subpop_models_ordered

    for field in fields(target):

        name = field.name

        # FluTravelStateTensors has an attribute
        #   that is a dictionary called `init_vals` --
        #   disregard, as this only used to store
        #   initial values for resetting, but is not
        #   used in the travel model computation
        if name == "init_vals":
            continue

        metapop_vals = []

        for model in subpop_models_ordered.values():
            current_val = getattr(model.state, name)
            metapop_vals.append(current_val)

        # Probably want to update this to be cleaner...
        # `SubpopState` fields that correspond to `Schedule` instances
        # have initial values of `None` -- but we cannot build a tensor
        # with `None` values, so we convert values to 0s.
        if any(v is None for v in metapop_vals):
            setattr(target, name, torch.tensor(np.full(np.shape(metapop_vals), 0.0)))
        else:
            setattr(target, name, torch.tensor(np.asarray(metapop_vals)))

update_travel_params_tensors() -> None

Update travel_params_tensors attribute in-place. FluTravelParamsTensors only has fields corresponding to parameters relevant for the travel model. Converts subpopulation-specific parameters to tensors of size (L, A, R) for location-age-risk (except for a few exceptions that have different dimensions).

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_travel_params_tensors(self) -> None:
    """
    Update `travel_params_tensors` attribute in-place.
    `FluTravelParamsTensors` only has fields corresponding
    to parameters relevant for the travel model.
    Converts subpopulation-specific parameters to
    tensors of size (L, A, R) for location-age-risk
    (except for a few exceptions that have different dimensions).
    """

    self.update_params_tensors(target=self.travel_params_tensors)

update_travel_state_tensors() -> None

Update travel_state_tensors attribute in-place. FluTravelStateTensors only has fields corresponding to state variables relevant for the travel model. Converts subpopulation-specific state to tensors of size (L, A, R) for location-age-risk (except for a few exceptions that have different dimensions).

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_travel_state_tensors(self) -> None:
    """
    Update `travel_state_tensors` attribute in-place.
    `FluTravelStateTensors` only has fields corresponding
    to state variables relevant for the travel model.
    Converts subpopulation-specific state to
    tensors of size (L, A, R) for location-age-risk
    (except for a few exceptions that have different dimensions).
    """

    self.update_state_tensors(self.travel_state_tensors)

FluMetapopModelError

Bases: MetapopModelError

Custom exceptions for flu metapopulation simulation model errors.

Source code in CLT_BaseModel/flu_core/flu_components.py
class FluMetapopModelError(clt.MetapopModelError):
    """Custom exceptions for flu metapopulation simulation model errors."""
    pass

FluMixingParams dataclass

Contains parameters corresponding to inter-subpopulation (metapopulation model) specifications: the number of subpopulations included, and the travel proportions between them.

Parameters:

Name Type Description Default
num_locations int

Number of locations (subpopulations) in the metapopulation model.

required
travel_proportions np.ndarray of shape (A, R)

L x L array of floats in [0,1], where L is the number of locations (subpopulations), and the i-jth element is the proportion of people in subpopulation i that travel to subpopulation j.

required
Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass(frozen=True)
class FluMixingParams:
    """
    Contains parameters corresponding to inter-subpopulation
    (metapopulation model) specifications: the number of
    subpopulations included, and the travel proportions between them.

    Params:
        num_locations (int):
            Number of locations (subpopulations) in the
            metapopulation model.
        travel_proportions (np.ndarray of shape (A, R)):
            L x L array of floats in [0,1], where L is the number
            of locations (subpopulations), and the i-jth element
            is the proportion of people in subpopulation i that
            travel to subpopulation j.
    """

    num_locations: Optional[int]
    travel_proportions: Optional[np.ndarray]

FluPrecomputedTensors

Stores precomputed quantities that are repeatedly used, for computational efficiency.

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
class FluPrecomputedTensors:
    """
    Stores precomputed quantities that are repeatedly
    used, for computational efficiency.
    """

    def __init__(self,
                 total_pop_LAR_tensor: torch.Tensor,
                 params: FluTravelParamsTensors) -> None:

        self.total_pop_LAR_tensor = total_pop_LAR_tensor

        self.L = int(params.num_locations.item())
        self.A = int(params.num_age_groups.item())
        self.R = int(params.num_risk_groups.item())

        self.total_pop_LA = torch.sum(self.total_pop_LAR_tensor, dim=2)

        # Remove the diagonal!
        self.nonlocal_travel_prop = params.travel_proportions.clone().fill_diagonal_(0.0)

        # We don't need einsum for residents traveling
        #   -- Dave and Remy helped me check this
        # \sum_{k \not = \ell} v^{\ell \rightarrow k}
        # Note we already have k \not = \ell because we set the diagonal of
        #   nonlocal_travel_prop to 0
        self.sum_residents_nonlocal_travel_prop = self.nonlocal_travel_prop.sum(dim=1)

FluSubpopModel

Bases: SubpopModel

Class for creating ImmunoSEIRS flu model with predetermined fixed structure -- initial values and epidemiological structure are populated by user-specified JSON files.

Key method create_transmission_model returns a SubpopModel instance with S-E-I-H-R-D compartments and M and MV epi metrics.

The update structure is as follows
  • S <- S + R_to_S - S_to_E
  • E <- E + S_to_E - E_to_IP - E_to_IA
  • IA <- IA + E_to_IA - IA_to_R
  • IP <- IP + E_to_IP - IP_to_ISR - IP_to_ISH
  • ISR <- ISR + IP_to_ISR - ISR_to_R
  • ISH <- ISH + IP_to_ISH - ISH_to_HR - ISH_to_HD
  • HR <- HR + ISH_to_HR - HR_to_R
  • HD <- HD + ISH_to_HD - HD_to_D
  • R <- R + ISR_to_R + HR_to_R - R_to_S
  • D <- D + HD_to_D
The following are TransitionVariable instances
  • R_to_S is a RecoveredToSusceptible instance
  • S_to_E is a SusceptibleToExposed instance
  • IP_to_ISR is a PresympToSympRecover instance
  • IP_to_ISH is a PresympToSympHospital instance
  • ISH_to_HR is a SympHospitalToHospRecover instance
  • ISH_to_HD is a SympHospitalToHospDead instance
  • ISR_to_R is a SympRecoverToRecovered instance
  • HR_to_R is a HospRecoverToRecovered instance
  • HD_to_D is a HospDeadToDead instance
There are three TransitionVariableGroups
  • E_out (handles E_to_IP and E_to_IA)
  • IP_out (handles IP_to_ISR and IP_to_ISH)
  • ISH_out (handles ISH_to_HR and ISH_to_HD)
The following are EpiMetric instances
  • M is a InfInducedImmunity instance
  • MV is a VaxInducedImmunity instance

Transition rates and update formulas are specified in corresponding classes.

See parent class SubpopModel's docstring for additional attributes.

Source code in CLT_BaseModel/flu_core/flu_components.py
1069
1070
1071
1072
1073
1074
1075
1076
1077
1078
1079
1080
1081
1082
1083
1084
1085
1086
1087
1088
1089
1090
1091
1092
1093
1094
1095
1096
1097
1098
1099
1100
1101
1102
1103
1104
1105
1106
1107
1108
1109
1110
1111
1112
1113
1114
1115
1116
1117
1118
1119
1120
1121
1122
1123
1124
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170
1171
1172
1173
1174
1175
1176
1177
1178
1179
1180
1181
1182
1183
1184
1185
1186
1187
1188
1189
1190
1191
1192
1193
1194
1195
1196
1197
1198
1199
1200
1201
1202
1203
1204
1205
1206
1207
1208
1209
1210
1211
1212
1213
1214
1215
1216
1217
1218
1219
1220
1221
1222
1223
1224
1225
1226
1227
1228
1229
1230
1231
1232
1233
1234
1235
1236
1237
1238
1239
1240
1241
1242
1243
1244
1245
1246
1247
1248
1249
1250
1251
1252
1253
1254
1255
1256
1257
1258
1259
1260
1261
1262
1263
1264
1265
1266
1267
1268
1269
1270
1271
1272
1273
1274
1275
1276
1277
1278
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289
1290
1291
1292
1293
1294
1295
1296
1297
1298
1299
1300
1301
1302
1303
1304
1305
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325
1326
1327
1328
1329
1330
1331
1332
1333
1334
1335
1336
1337
1338
1339
1340
1341
1342
1343
1344
1345
1346
1347
1348
1349
1350
1351
1352
1353
1354
1355
1356
1357
1358
1359
1360
1361
1362
1363
1364
1365
1366
1367
1368
1369
1370
1371
1372
1373
1374
1375
1376
1377
1378
1379
1380
1381
1382
1383
1384
1385
1386
1387
1388
1389
1390
1391
1392
1393
1394
1395
1396
1397
1398
1399
1400
1401
1402
1403
1404
1405
1406
1407
1408
1409
1410
1411
1412
1413
1414
1415
1416
1417
1418
1419
1420
1421
1422
1423
1424
1425
1426
1427
1428
1429
1430
1431
1432
1433
1434
1435
1436
1437
1438
1439
1440
1441
1442
1443
1444
1445
1446
1447
1448
1449
1450
1451
1452
1453
1454
1455
1456
1457
1458
1459
1460
1461
1462
1463
1464
1465
1466
1467
1468
1469
1470
1471
1472
1473
1474
1475
1476
1477
1478
1479
1480
1481
1482
1483
1484
1485
1486
1487
1488
1489
1490
1491
1492
1493
1494
1495
1496
1497
1498
1499
1500
1501
1502
1503
1504
1505
1506
1507
1508
1509
1510
1511
1512
1513
1514
1515
1516
1517
1518
1519
1520
1521
1522
1523
1524
1525
1526
1527
1528
1529
1530
1531
1532
1533
1534
1535
1536
1537
1538
1539
1540
1541
1542
1543
1544
1545
1546
1547
1548
1549
1550
1551
1552
1553
1554
1555
1556
1557
1558
1559
1560
1561
1562
1563
1564
1565
1566
1567
1568
1569
1570
1571
1572
1573
1574
1575
1576
1577
1578
1579
1580
1581
1582
1583
1584
1585
1586
class FluSubpopModel(clt.SubpopModel):
    """
    Class for creating ImmunoSEIRS flu model with predetermined fixed
    structure -- initial values and epidemiological structure are
    populated by user-specified `JSON` files.

    Key method create_transmission_model returns a `SubpopModel`
    instance with S-E-I-H-R-D compartments and M
    and MV epi metrics.

    The update structure is as follows:
        - S <- S + R_to_S - S_to_E
        - E <- E + S_to_E - E_to_IP - E_to_IA
        - IA <- IA + E_to_IA - IA_to_R 
        - IP <- IP + E_to_IP - IP_to_ISR - IP_to_ISH
        - ISR <- ISR + IP_to_ISR - ISR_to_R
        - ISH <- ISH + IP_to_ISH - ISH_to_HR - ISH_to_HD
        - HR <- HR + ISH_to_HR - HR_to_R
        - HD <- HD + ISH_to_HD - HD_to_D
        - R <- R + ISR_to_R + HR_to_R - R_to_S
        - D <- D + HD_to_D

    The following are TransitionVariable instances:
        - R_to_S is a RecoveredToSusceptible instance
        - S_to_E is a SusceptibleToExposed instance
        - IP_to_ISR is a PresympToSympRecover instance
        - IP_to_ISH is a PresympToSympHospital instance
        - ISH_to_HR is a SympHospitalToHospRecover instance
        - ISH_to_HD is a SympHospitalToHospDead instance
        - ISR_to_R is a SympRecoverToRecovered instance
        - HR_to_R is a HospRecoverToRecovered instance 
        - HD_to_D is a HospDeadToDead instance

    There are three TransitionVariableGroups:
        - E_out (handles E_to_IP and E_to_IA)
        - IP_out (handles IP_to_ISR and IP_to_ISH)
        - ISH_out (handles ISH_to_HR and ISH_to_HD)

    The following are EpiMetric instances:
        - M is a InfInducedImmunity instance
        - MV is a VaxInducedImmunity instance

    Transition rates and update formulas are specified in
    corresponding classes.

    See parent class `SubpopModel`'s docstring for additional attributes.
    """

    def __init__(self,
                 state: FluSubpopState,
                 params: FluSubpopParams,
                 simulation_settings: FluSubpopSchedules,
                 RNG: np.random.Generator,
                 schedules_spec: FluSubpopSchedules,
                 name: str):
        """
        Args:
            state (FluSubpopState):
                holds current simulation state information,
                such as current values of epidemiological compartments
                and epi metrics.
            params (FluSubpopParams):
                holds epidemiological parameter values.
            simulation_settings (SimulationSettings):
                holds simulation settings.
            RNG (np.random.Generator):
                numpy random generator object used to obtain
                random numbers.
            schedules_spec (FluSubpopSchedules):
                holds dataframes that specify `Schedule` instances.
            name (str):
                unique name of MetapopModel instance.
        """

        self.schedules_spec = schedules_spec

        # IMPORTANT NOTE: as always, we must be careful with mutable objects
        # and generally use deep copies to avoid modification of the same
        # object. But in this function call, using deep copies is unnecessary
        # (redundant) because the parent class `SubpopModel`'s `__init__`
        # creates deep copies.
        super().__init__(state, params, simulation_settings, RNG, name)

        self.params = clt.updated_dataclass(self.params, {"start_real_date": self.start_real_date})

        self.update_vax_induced_risk_reduce_initial()

    def update_vax_induced_risk_reduce_initial(self) -> None:
        """
        Recomputes `vax_induced_inf_risk_reduce_initial`,
        `vax_induced_hosp_risk_reduce_initial`, and
        `vax_induced_death_risk_reduce_initial` from the current
        `daily_vaccines` schedule and updates `self.params` in place.

        This must be re-run (not just computed once at construction)
        whenever the `daily_vaccines` schedule or any of the underlying
        base parameters (`vax_induced_*_risk_reduce`,
        `vax_induced_immune_wane`) change after construction -- e.g. via
        `replace_schedule` or a `ScenarioRunner` parameter override --
        otherwise these derived values would silently keep reflecting
        the schedule/params from construction time. See
        `reset_simulation`, which calls this for the same reason
        `VaxInducedImmunity`'s initial value is recomputed there.
        """

        inf_initial, hosp_initial, death_initial = compute_vax_induced_risk_reduce_initial(
            self.params, self.schedules, self.start_real_date)
        self.params = clt.updated_dataclass(self.params, {
            "vax_induced_inf_risk_reduce_initial": inf_initial,
            "vax_induced_hosp_risk_reduce_initial": hosp_initial,
            "vax_induced_death_risk_reduce_initial": death_initial,
        })

    def check_humidity_input(self) -> None:
        """
        Check that absolute humidity values are non-negative.
        """

        humidity_values = self.schedules['absolute_humidity'].timeseries_df['absolute_humidity'].values
        if np.any(humidity_values < 0):
            raise FluSubpopModelError("Error: absolute humidity values must be non-negative.")

    def check_vaccination_input(self) -> None:
        """
        Check that vaccination values are positive.
        If vaccinations exceed 100% over a year, issue a warning.
        """

        df_vaccine = self.schedules['daily_vaccines'].timeseries_df.copy()

        ## Check all entries are positive
        all_positive = all([
            (x >= 0).all() for x in df_vaccine['daily_vaccines'].values
            ])
        if not(all_positive):
            raise FluSubpopModelError("Error: vaccination values must be non-negative.")

        ## Check cumulative vaccination never exceeds 100% over 365 days
        df_vaccine['datetime'] = pd.to_datetime(df_vaccine.index)
        df_vaccine.set_index('datetime', inplace=True)

        # Ensure there is one row per day within time range
        full_date_range = pd.date_range(start=df_vaccine.index.min(), end=df_vaccine.index.max(), freq='D')
        df_vaccine = df_vaccine.reindex(full_date_range)

        # Replace NAs due to reindexing (if any)
        na_rows = df_vaccine['daily_vaccines'].isna()
        value_shape = df_vaccine['daily_vaccines'].values[0].shape
        df_vaccine.loc[na_rows, 'daily_vaccines'] = \
            pd.Series(
                [np.zeros(value_shape)] * na_rows.sum(), 
                index=df_vaccine.loc[na_rows].index
                )

        # Calculate rolling 1-year sum of vaccination rates
        window_size_days = min(365, len(df_vaccine))

        data_windows = np.lib.stride_tricks.sliding_window_view(
            df_vaccine['daily_vaccines'].values, 
            window_size_days
            )
        vaccines_rolling_sum = np.sum(data_windows, axis=-1)

        # Check whether any cumulative vaccinations exceed 100%
        max_values_above_one = [x.max() > 1 for x in vaccines_rolling_sum]

        # Find first index rolling sum exceeds 100% for some age group (if there is one)
        if sum(max_values_above_one) > 0:
            first_exceeds_idx = max_values_above_one.index(True)

            vaccines_cml_exceeds = vaccines_rolling_sum[first_exceeds_idx]
            exceeds_first_date = df_vaccine.index[first_exceeds_idx + window_size_days - 1]

            msg = 'Cumulative vaccination over a 365-day period exceeds 100% on (at least) ' +\
                f'the following date: {exceeds_first_date}. Cumulative vaccination by that date is \n' +\
                str(vaccines_cml_exceeds)
            warnings.warn(msg)

    def check_calendar_variables_input(self) -> None:
        """
        Check school and calendar variables in flu_contact_matrix
        schedule are between 0 and 1
        """

        flu_contact_matrix = self.schedules['flu_contact_matrix'].timeseries_df

        for variable in ['is_school_day', 'is_work_day']:
            values = flu_contact_matrix[variable].values

            if (values.min() < 0) or (values.max() > 1):
                msg = f'Error: {variable} values must be between 0 and 1.'
                raise FluSubpopModelError(msg)

    def check_contact_matrix_input(self) -> None:
        """
        Check contact matrix entries are non-negative.
        Check total contact is greater than the sum of the
        school and work matrices.
        """   

        if not(np.all(self.params.total_contact_matrix >= 0)):
            raise FluSubpopModelError(
                'Some entries of the total contact matrix are negative \n'+\
                f'{self.params.total_contact_matrix}'
                )

        if not(np.all(self.params.school_contact_matrix >= 0)):
            raise FluSubpopModelError(
                'Some entries of the school contact matrix are negative \n'+\
                f'{self.params.school_contact_matrix}'
                )

        if not(np.all(self.params.work_contact_matrix >= 0)):
            raise FluSubpopModelError(
                'Some entries of the work contact matrix are negative \n'+\
                f'{self.params.work_contact_matrix}'
                )

        if not(np.all((
            self.params.total_contact_matrix - self.params.school_contact_matrix - 
            self.params.work_contact_matrix) >= 0)):
            raise FluSubpopModelError(
                'The total contact matrix must be at least greater than the sum of ' +\
                'the work and school contact matrices.'
            )

    def check_rate_input(self) -> None:
        """
        Ensure all rate values are strictly positive, and other
        variables (waning, saturation, reductions) are non-negative.
        """

        p = self.params
        rates_list = [
            p.R_to_S_rate, p.E_to_I_rate, p.IP_to_IS_rate, p.ISH_to_H_rate,
            p.ISR_to_R_rate, p.IA_to_R_rate, p.HR_to_R_rate, p.HD_to_D_rate,
            p.E_to_IA_prop]

        other_params_list = [
            p.humidity_impact, p.inf_induced_saturation, p.inf_induced_immune_wane,
            p.vax_induced_saturation, p.vax_induced_immune_wane, p.inf_induced_inf_risk_reduce,
            p.inf_induced_hosp_risk_reduce, p.inf_induced_death_risk_reduce, 
            p.vax_induced_inf_risk_reduce, p.vax_induced_hosp_risk_reduce,
            p.vax_induced_death_risk_reduce, p.IP_relative_inf,
            p.IA_relative_inf, p.relative_suscept,
            p.ISH_to_HD_prop, p.IP_to_ISH_prop, p.beta_baseline
        ]

        for value in rates_list:
            if not(np.all(value >= 0)):
                raise FluSubpopModelError('All transition rates must be positive values.')
            if not(np.all(value > 0)):
                msg = 'Some transition rates are equal to zero.'
                warnings.warn(msg)

        for value in other_params_list:
            if not(np.all(value >= 0)):
                raise FluSubpopModelError('Some parameter values are negative.')

    def check_initial_compartment_input(self) -> None:
        """
        Ensure all initial compartment and saturation values are non-negative.
        """

        compartments_list = ["S", "E", "IP", "ISR", "ISH", "IA", "HR", "HD", "R", "D"]
        immunity_values_list = ["M", "MV"]

        for state_name in compartments_list + immunity_values_list:
            if not(np.all(getattr(self.state, state_name) >= 0)):
                raise FluSubpopModelError(
                    'Initial compartment and immunity values must be non-negative. ' +\
                    f'{state_name} is negative: {getattr(self.state, state_name)} for subpopulation ' +\
                    f'{self.name}.'
                )

    def run_input_checks(self) -> None:
        """
        Check the following:
            - if total vaccinations exceed 100% over a year we issue a warning
            - school and work calendar variables must be between 0 and 1
            - absolute humidity values, contact matrix entries, daily vaccination
              must be non-negative
            - total contact matrix must be greater than the sum of the school
              and work contact matrices
            - all rate values must be strictly positive
            - initial compartmental values must be non-negative
        """

        self.check_humidity_input()
        self.check_vaccination_input()
        self.check_calendar_variables_input()
        self.check_contact_matrix_input()
        self.check_rate_input()
        self.check_initial_compartment_input()

    def prepare_daily_state(self) -> None:
        """
        Override parent method to add vaccine immunity reset check.
        At beginning of each day, update schedules, dynamic values,
        and check for vaccine immunity reset.
        """
        # Call parent implementation first to update schedules and dynamic vals
        super().prepare_daily_state()

        # Check and potentially reset vaccine-induced immunity
        if hasattr(self.epi_metrics, 'MV'):
            self.epi_metrics.MV.check_and_apply_reset(
                self.current_real_date,
                self.params
            )

    def create_compartments(self) -> sc.objdict[str, clt.Compartment]:

        # Create `Compartment` instances S-E-IA-IP-IS-H-R-D (7 compartments total)
        # Save instances in `sc.objdict` and return objdict

        compartments = sc.objdict()

        for name in ("S", "E", "IP", "ISR", "ISH", "IA", "HR", "HD", "R", "D"):
            compartments[name] = clt.Compartment(getattr(self.state, name))

        return compartments

    def create_dynamic_vals(self) -> sc.objdict[str, clt.DynamicVal]:
        """
        Create all `DynamicVal` instances, save in `sc.objdict`, and return objdict
        """

        dynamic_vals = sc.objdict()

        dynamic_vals["beta_reduce"] = BetaReduce(init_val=0.0,
                                                 is_enabled=False)

        return dynamic_vals

    def create_schedules(self) -> sc.objdict[str, clt.Schedule]:
        """
        Create all `Schedule` instances, save in `sc.objdict`, and return objdict
        """

        schedules = sc.objdict()

        schedules["absolute_humidity"] = AbsoluteHumidity()
        schedules["flu_contact_matrix"] = FluContactMatrix()
        schedules["daily_vaccines"] = DailyVaccines(
            vax_protection_delay_days=self.params.vax_protection_delay_days
        )
        schedules["mobility_modifier"] = MobilityModifier()

        for field, df in asdict(self.schedules_spec).items():

            try:
                if 'day_of_week' not in df.columns:
                    df["date"] = pd.to_datetime(df["date"], format='%Y-%m-%d').dt.date
            except ValueError as e:
                raise ValueError("Error: dates should be strings in YYYY-MM-DD format or "
                                 "`date.datetime` objects.") from e

            schedules[field].timeseries_df = df
            schedules[field].postprocess_data_input()

        return schedules

    def create_transition_variables(self) -> sc.objdict[str, clt.TransitionVariable]:
        """
        Create all `TransitionVariable` instances,
        save in `sc.objdict`, and return objdict
        """

        # NOTE: see the parent class `SubpopModel`'s `__init__` --
        # `create_transition_variables` is called after
        # `simulation_settings` is assigned

        transition_type = self.simulation_settings.transition_type

        transition_variables = sc.objdict()

        S = self.compartments.S
        E = self.compartments.E
        IP = self.compartments.IP
        ISR = self.compartments.ISR
        ISH = self.compartments.ISH
        IA = self.compartments.IA
        HR = self.compartments.HR
        HD = self.compartments.HD
        R = self.compartments.R
        D = self.compartments.D

        transition_variables.R_to_S = RecoveredToSusceptible(R, S, transition_type)
        transition_variables.S_to_E = SusceptibleToExposed(S, E, transition_type)
        transition_variables.IP_to_ISR = PresympToSympRecover(IP, ISR, transition_type, True)
        transition_variables.IP_to_ISH = PresympToSympHospital(IP, ISH, transition_type, True)
        transition_variables.IA_to_R = AsympToRecovered(IA, R, transition_type)
        transition_variables.E_to_IP = ExposedToPresymp(E, IP, transition_type, True)
        transition_variables.E_to_IA = ExposedToAsymp(E, IA, transition_type, True)
        transition_variables.ISR_to_R = SympRecoverToRecovered(ISR, R, transition_type)
        transition_variables.ISH_to_HR = SympHospitalToHospRecover(ISH, HR, transition_type, True)
        transition_variables.ISH_to_HD = SympHospitalToHospDead(ISH, HD, transition_type, True)
        transition_variables.HR_to_R = HospRecoverToRecovered(HR, R, transition_type)
        transition_variables.HD_to_D = HospDeadToDead(HD, D, transition_type)

        return transition_variables

    def create_transition_variable_groups(self) -> sc.objdict[str, clt.TransitionVariableGroup]:
        """
        Create all transition variable groups described in docstring (3 transition
        variable groups total), save in `sc.objdict`, return objdict
        """

        # Shortcuts for attribute access
        # NOTE: see the parent class `SubpopModel`'s `__init__` --
        # `create_transition_variable_groups` is called after
        # `simulation_settings` is assigned

        transition_type = self.simulation_settings.transition_type

        transition_variable_groups = sc.objdict()

        transition_variable_groups.E_out = clt.TransitionVariableGroup(self.compartments.E,
                                                                       transition_type,
                                                                       (self.transition_variables.E_to_IP,
                                                                        self.transition_variables.E_to_IA))

        transition_variable_groups.IP_out = clt.TransitionVariableGroup(self.compartments.IP,
                                                                        transition_type,
                                                                        (self.transition_variables.IP_to_ISR,
                                                                         self.transition_variables.IP_to_ISH))

        transition_variable_groups.ISH_out = clt.TransitionVariableGroup(self.compartments.ISH,
                                                                         transition_type,
                                                                         (self.transition_variables.ISH_to_HR,
                                                                          self.transition_variables.ISH_to_HD))

        return transition_variable_groups

    def create_epi_metrics(self) -> sc.objdict[str, clt.EpiMetric]:
        """
        Create all epi metric described in docstring (2 state
        variables total), save in `sc.objdict`, and return objdict
        """

        epi_metrics = sc.objdict()

        epi_metrics.M = \
            InfInducedImmunity(getattr(self.state, "M"),
                               self.transition_variables.R_to_S)

        epi_metrics.MV = \
            VaxInducedImmunity(getattr(self.state, "MV"),
                               self.current_real_date,
                               self.params,
                               self.schedules,
                               self.simulation_settings.timesteps_per_day)

        return epi_metrics

    def modify_subpop_params(self,
                             updates_dict: dict):
        """
        This method lets users safely modify a single subpopulation
        parameters field; if this subpop model is associated with
        a metapop model, the metapopulation-wide tensors are updated
        automatically afterward. See also `modify_subpop_params` method on
        `FluMetapopModel`.

        Parameters:
            updates_dict (dict):
                Dictionary specifying values to update in a
                `FluSubpopParams` instance -- keys must match the
                field names of `FluSubpopParams`.
        """

        # If associated with metapop model, run this method
        #   on the metapop model itself to handle metapopulation-wide
        #   tensor updating
        if self.metapop_model:
            self.metapop_model.modify_subpop_params(self.name,
                                                    updates_dict)
        else:
            # Since `SubpopParams` is frozen, we return a new instance
            #   with the reflected updates
            self.params = clt.updated_dataclass(self.params, updates_dict)

    def reset_simulation(self) -> None:
        """
        Extends the base `reset_simulation` to recompute `MV.init_val` and
        `vax_induced_*_risk_reduce_initial` from the currently loaded
        vaccine schedule (and current base params) before resetting.

        This ensures that if the `daily_vaccines` schedule has been replaced
        (e.g. via `replace_schedule`), or the underlying
        `vax_induced_*_risk_reduce`/`vax_induced_immune_wane` params have
        been overridden (e.g. by `ScenarioRunner`), the model resets to
        values consistent with the current schedule/params, rather than
        the values computed at construction time from the original
        schedule/params.

        The `MV.init_val` recomputation uses
        `VaxInducedImmunity.adjust_initial_value()` with
        `MV.original_init_val` as the base β€” the unmodified value
        from the state JSON β€” so adjustments do not compound across calls.
        """

        MV = self.epi_metrics["MV"]
        new_init_val = MV.adjust_initial_value(
            MV.original_init_val,
            self.start_real_date,
            self.params,
            self.schedules,
            self.simulation_settings.timesteps_per_day,
        )
        # Use the init_val setter so current_val is also updated immediately,
        # before super()'s reset loop overwrites it again (harmlessly).
        MV.init_val = new_init_val

        self.update_vax_induced_risk_reduce_initial()

        super().reset_simulation()

__init__(state: FluSubpopState, params: FluSubpopParams, simulation_settings: FluSubpopSchedules, RNG: np.random.Generator, schedules_spec: FluSubpopSchedules, name: str)

Parameters:

Name Type Description Default
state FluSubpopState

holds current simulation state information, such as current values of epidemiological compartments and epi metrics.

required
params FluSubpopParams

holds epidemiological parameter values.

required
simulation_settings SimulationSettings

holds simulation settings.

required
RNG Generator

numpy random generator object used to obtain random numbers.

required
schedules_spec FluSubpopSchedules

holds dataframes that specify Schedule instances.

required
name str

unique name of MetapopModel instance.

required
Source code in CLT_BaseModel/flu_core/flu_components.py
def __init__(self,
             state: FluSubpopState,
             params: FluSubpopParams,
             simulation_settings: FluSubpopSchedules,
             RNG: np.random.Generator,
             schedules_spec: FluSubpopSchedules,
             name: str):
    """
    Args:
        state (FluSubpopState):
            holds current simulation state information,
            such as current values of epidemiological compartments
            and epi metrics.
        params (FluSubpopParams):
            holds epidemiological parameter values.
        simulation_settings (SimulationSettings):
            holds simulation settings.
        RNG (np.random.Generator):
            numpy random generator object used to obtain
            random numbers.
        schedules_spec (FluSubpopSchedules):
            holds dataframes that specify `Schedule` instances.
        name (str):
            unique name of MetapopModel instance.
    """

    self.schedules_spec = schedules_spec

    # IMPORTANT NOTE: as always, we must be careful with mutable objects
    # and generally use deep copies to avoid modification of the same
    # object. But in this function call, using deep copies is unnecessary
    # (redundant) because the parent class `SubpopModel`'s `__init__`
    # creates deep copies.
    super().__init__(state, params, simulation_settings, RNG, name)

    self.params = clt.updated_dataclass(self.params, {"start_real_date": self.start_real_date})

    self.update_vax_induced_risk_reduce_initial()

check_calendar_variables_input() -> None

Check school and calendar variables in flu_contact_matrix schedule are between 0 and 1

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_calendar_variables_input(self) -> None:
    """
    Check school and calendar variables in flu_contact_matrix
    schedule are between 0 and 1
    """

    flu_contact_matrix = self.schedules['flu_contact_matrix'].timeseries_df

    for variable in ['is_school_day', 'is_work_day']:
        values = flu_contact_matrix[variable].values

        if (values.min() < 0) or (values.max() > 1):
            msg = f'Error: {variable} values must be between 0 and 1.'
            raise FluSubpopModelError(msg)

check_contact_matrix_input() -> None

Check contact matrix entries are non-negative. Check total contact is greater than the sum of the school and work matrices.

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_contact_matrix_input(self) -> None:
    """
    Check contact matrix entries are non-negative.
    Check total contact is greater than the sum of the
    school and work matrices.
    """   

    if not(np.all(self.params.total_contact_matrix >= 0)):
        raise FluSubpopModelError(
            'Some entries of the total contact matrix are negative \n'+\
            f'{self.params.total_contact_matrix}'
            )

    if not(np.all(self.params.school_contact_matrix >= 0)):
        raise FluSubpopModelError(
            'Some entries of the school contact matrix are negative \n'+\
            f'{self.params.school_contact_matrix}'
            )

    if not(np.all(self.params.work_contact_matrix >= 0)):
        raise FluSubpopModelError(
            'Some entries of the work contact matrix are negative \n'+\
            f'{self.params.work_contact_matrix}'
            )

    if not(np.all((
        self.params.total_contact_matrix - self.params.school_contact_matrix - 
        self.params.work_contact_matrix) >= 0)):
        raise FluSubpopModelError(
            'The total contact matrix must be at least greater than the sum of ' +\
            'the work and school contact matrices.'
        )

check_humidity_input() -> None

Check that absolute humidity values are non-negative.

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_humidity_input(self) -> None:
    """
    Check that absolute humidity values are non-negative.
    """

    humidity_values = self.schedules['absolute_humidity'].timeseries_df['absolute_humidity'].values
    if np.any(humidity_values < 0):
        raise FluSubpopModelError("Error: absolute humidity values must be non-negative.")

check_initial_compartment_input() -> None

Ensure all initial compartment and saturation values are non-negative.

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_initial_compartment_input(self) -> None:
    """
    Ensure all initial compartment and saturation values are non-negative.
    """

    compartments_list = ["S", "E", "IP", "ISR", "ISH", "IA", "HR", "HD", "R", "D"]
    immunity_values_list = ["M", "MV"]

    for state_name in compartments_list + immunity_values_list:
        if not(np.all(getattr(self.state, state_name) >= 0)):
            raise FluSubpopModelError(
                'Initial compartment and immunity values must be non-negative. ' +\
                f'{state_name} is negative: {getattr(self.state, state_name)} for subpopulation ' +\
                f'{self.name}.'
            )

check_rate_input() -> None

Ensure all rate values are strictly positive, and other variables (waning, saturation, reductions) are non-negative.

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_rate_input(self) -> None:
    """
    Ensure all rate values are strictly positive, and other
    variables (waning, saturation, reductions) are non-negative.
    """

    p = self.params
    rates_list = [
        p.R_to_S_rate, p.E_to_I_rate, p.IP_to_IS_rate, p.ISH_to_H_rate,
        p.ISR_to_R_rate, p.IA_to_R_rate, p.HR_to_R_rate, p.HD_to_D_rate,
        p.E_to_IA_prop]

    other_params_list = [
        p.humidity_impact, p.inf_induced_saturation, p.inf_induced_immune_wane,
        p.vax_induced_saturation, p.vax_induced_immune_wane, p.inf_induced_inf_risk_reduce,
        p.inf_induced_hosp_risk_reduce, p.inf_induced_death_risk_reduce, 
        p.vax_induced_inf_risk_reduce, p.vax_induced_hosp_risk_reduce,
        p.vax_induced_death_risk_reduce, p.IP_relative_inf,
        p.IA_relative_inf, p.relative_suscept,
        p.ISH_to_HD_prop, p.IP_to_ISH_prop, p.beta_baseline
    ]

    for value in rates_list:
        if not(np.all(value >= 0)):
            raise FluSubpopModelError('All transition rates must be positive values.')
        if not(np.all(value > 0)):
            msg = 'Some transition rates are equal to zero.'
            warnings.warn(msg)

    for value in other_params_list:
        if not(np.all(value >= 0)):
            raise FluSubpopModelError('Some parameter values are negative.')

check_vaccination_input() -> None

Check that vaccination values are positive. If vaccinations exceed 100% over a year, issue a warning.

Source code in CLT_BaseModel/flu_core/flu_components.py
def check_vaccination_input(self) -> None:
    """
    Check that vaccination values are positive.
    If vaccinations exceed 100% over a year, issue a warning.
    """

    df_vaccine = self.schedules['daily_vaccines'].timeseries_df.copy()

    ## Check all entries are positive
    all_positive = all([
        (x >= 0).all() for x in df_vaccine['daily_vaccines'].values
        ])
    if not(all_positive):
        raise FluSubpopModelError("Error: vaccination values must be non-negative.")

    ## Check cumulative vaccination never exceeds 100% over 365 days
    df_vaccine['datetime'] = pd.to_datetime(df_vaccine.index)
    df_vaccine.set_index('datetime', inplace=True)

    # Ensure there is one row per day within time range
    full_date_range = pd.date_range(start=df_vaccine.index.min(), end=df_vaccine.index.max(), freq='D')
    df_vaccine = df_vaccine.reindex(full_date_range)

    # Replace NAs due to reindexing (if any)
    na_rows = df_vaccine['daily_vaccines'].isna()
    value_shape = df_vaccine['daily_vaccines'].values[0].shape
    df_vaccine.loc[na_rows, 'daily_vaccines'] = \
        pd.Series(
            [np.zeros(value_shape)] * na_rows.sum(), 
            index=df_vaccine.loc[na_rows].index
            )

    # Calculate rolling 1-year sum of vaccination rates
    window_size_days = min(365, len(df_vaccine))

    data_windows = np.lib.stride_tricks.sliding_window_view(
        df_vaccine['daily_vaccines'].values, 
        window_size_days
        )
    vaccines_rolling_sum = np.sum(data_windows, axis=-1)

    # Check whether any cumulative vaccinations exceed 100%
    max_values_above_one = [x.max() > 1 for x in vaccines_rolling_sum]

    # Find first index rolling sum exceeds 100% for some age group (if there is one)
    if sum(max_values_above_one) > 0:
        first_exceeds_idx = max_values_above_one.index(True)

        vaccines_cml_exceeds = vaccines_rolling_sum[first_exceeds_idx]
        exceeds_first_date = df_vaccine.index[first_exceeds_idx + window_size_days - 1]

        msg = 'Cumulative vaccination over a 365-day period exceeds 100% on (at least) ' +\
            f'the following date: {exceeds_first_date}. Cumulative vaccination by that date is \n' +\
            str(vaccines_cml_exceeds)
        warnings.warn(msg)

create_dynamic_vals() -> sc.objdict[str, clt.DynamicVal]

Create all DynamicVal instances, save in sc.objdict, and return objdict

Source code in CLT_BaseModel/flu_core/flu_components.py
def create_dynamic_vals(self) -> sc.objdict[str, clt.DynamicVal]:
    """
    Create all `DynamicVal` instances, save in `sc.objdict`, and return objdict
    """

    dynamic_vals = sc.objdict()

    dynamic_vals["beta_reduce"] = BetaReduce(init_val=0.0,
                                             is_enabled=False)

    return dynamic_vals

create_epi_metrics() -> sc.objdict[str, clt.EpiMetric]

Create all epi metric described in docstring (2 state variables total), save in sc.objdict, and return objdict

Source code in CLT_BaseModel/flu_core/flu_components.py
def create_epi_metrics(self) -> sc.objdict[str, clt.EpiMetric]:
    """
    Create all epi metric described in docstring (2 state
    variables total), save in `sc.objdict`, and return objdict
    """

    epi_metrics = sc.objdict()

    epi_metrics.M = \
        InfInducedImmunity(getattr(self.state, "M"),
                           self.transition_variables.R_to_S)

    epi_metrics.MV = \
        VaxInducedImmunity(getattr(self.state, "MV"),
                           self.current_real_date,
                           self.params,
                           self.schedules,
                           self.simulation_settings.timesteps_per_day)

    return epi_metrics

create_schedules() -> sc.objdict[str, clt.Schedule]

Create all Schedule instances, save in sc.objdict, and return objdict

Source code in CLT_BaseModel/flu_core/flu_components.py
def create_schedules(self) -> sc.objdict[str, clt.Schedule]:
    """
    Create all `Schedule` instances, save in `sc.objdict`, and return objdict
    """

    schedules = sc.objdict()

    schedules["absolute_humidity"] = AbsoluteHumidity()
    schedules["flu_contact_matrix"] = FluContactMatrix()
    schedules["daily_vaccines"] = DailyVaccines(
        vax_protection_delay_days=self.params.vax_protection_delay_days
    )
    schedules["mobility_modifier"] = MobilityModifier()

    for field, df in asdict(self.schedules_spec).items():

        try:
            if 'day_of_week' not in df.columns:
                df["date"] = pd.to_datetime(df["date"], format='%Y-%m-%d').dt.date
        except ValueError as e:
            raise ValueError("Error: dates should be strings in YYYY-MM-DD format or "
                             "`date.datetime` objects.") from e

        schedules[field].timeseries_df = df
        schedules[field].postprocess_data_input()

    return schedules

create_transition_variable_groups() -> sc.objdict[str, clt.TransitionVariableGroup]

Create all transition variable groups described in docstring (3 transition variable groups total), save in sc.objdict, return objdict

Source code in CLT_BaseModel/flu_core/flu_components.py
def create_transition_variable_groups(self) -> sc.objdict[str, clt.TransitionVariableGroup]:
    """
    Create all transition variable groups described in docstring (3 transition
    variable groups total), save in `sc.objdict`, return objdict
    """

    # Shortcuts for attribute access
    # NOTE: see the parent class `SubpopModel`'s `__init__` --
    # `create_transition_variable_groups` is called after
    # `simulation_settings` is assigned

    transition_type = self.simulation_settings.transition_type

    transition_variable_groups = sc.objdict()

    transition_variable_groups.E_out = clt.TransitionVariableGroup(self.compartments.E,
                                                                   transition_type,
                                                                   (self.transition_variables.E_to_IP,
                                                                    self.transition_variables.E_to_IA))

    transition_variable_groups.IP_out = clt.TransitionVariableGroup(self.compartments.IP,
                                                                    transition_type,
                                                                    (self.transition_variables.IP_to_ISR,
                                                                     self.transition_variables.IP_to_ISH))

    transition_variable_groups.ISH_out = clt.TransitionVariableGroup(self.compartments.ISH,
                                                                     transition_type,
                                                                     (self.transition_variables.ISH_to_HR,
                                                                      self.transition_variables.ISH_to_HD))

    return transition_variable_groups

create_transition_variables() -> sc.objdict[str, clt.TransitionVariable]

Create all TransitionVariable instances, save in sc.objdict, and return objdict

Source code in CLT_BaseModel/flu_core/flu_components.py
def create_transition_variables(self) -> sc.objdict[str, clt.TransitionVariable]:
    """
    Create all `TransitionVariable` instances,
    save in `sc.objdict`, and return objdict
    """

    # NOTE: see the parent class `SubpopModel`'s `__init__` --
    # `create_transition_variables` is called after
    # `simulation_settings` is assigned

    transition_type = self.simulation_settings.transition_type

    transition_variables = sc.objdict()

    S = self.compartments.S
    E = self.compartments.E
    IP = self.compartments.IP
    ISR = self.compartments.ISR
    ISH = self.compartments.ISH
    IA = self.compartments.IA
    HR = self.compartments.HR
    HD = self.compartments.HD
    R = self.compartments.R
    D = self.compartments.D

    transition_variables.R_to_S = RecoveredToSusceptible(R, S, transition_type)
    transition_variables.S_to_E = SusceptibleToExposed(S, E, transition_type)
    transition_variables.IP_to_ISR = PresympToSympRecover(IP, ISR, transition_type, True)
    transition_variables.IP_to_ISH = PresympToSympHospital(IP, ISH, transition_type, True)
    transition_variables.IA_to_R = AsympToRecovered(IA, R, transition_type)
    transition_variables.E_to_IP = ExposedToPresymp(E, IP, transition_type, True)
    transition_variables.E_to_IA = ExposedToAsymp(E, IA, transition_type, True)
    transition_variables.ISR_to_R = SympRecoverToRecovered(ISR, R, transition_type)
    transition_variables.ISH_to_HR = SympHospitalToHospRecover(ISH, HR, transition_type, True)
    transition_variables.ISH_to_HD = SympHospitalToHospDead(ISH, HD, transition_type, True)
    transition_variables.HR_to_R = HospRecoverToRecovered(HR, R, transition_type)
    transition_variables.HD_to_D = HospDeadToDead(HD, D, transition_type)

    return transition_variables

modify_subpop_params(updates_dict: dict)

This method lets users safely modify a single subpopulation parameters field; if this subpop model is associated with a metapop model, the metapopulation-wide tensors are updated automatically afterward. See also modify_subpop_params method on FluMetapopModel.

Parameters:

Name Type Description Default
updates_dict dict

Dictionary specifying values to update in a FluSubpopParams instance -- keys must match the field names of FluSubpopParams.

required
Source code in CLT_BaseModel/flu_core/flu_components.py
def modify_subpop_params(self,
                         updates_dict: dict):
    """
    This method lets users safely modify a single subpopulation
    parameters field; if this subpop model is associated with
    a metapop model, the metapopulation-wide tensors are updated
    automatically afterward. See also `modify_subpop_params` method on
    `FluMetapopModel`.

    Parameters:
        updates_dict (dict):
            Dictionary specifying values to update in a
            `FluSubpopParams` instance -- keys must match the
            field names of `FluSubpopParams`.
    """

    # If associated with metapop model, run this method
    #   on the metapop model itself to handle metapopulation-wide
    #   tensor updating
    if self.metapop_model:
        self.metapop_model.modify_subpop_params(self.name,
                                                updates_dict)
    else:
        # Since `SubpopParams` is frozen, we return a new instance
        #   with the reflected updates
        self.params = clt.updated_dataclass(self.params, updates_dict)

prepare_daily_state() -> None

Override parent method to add vaccine immunity reset check. At beginning of each day, update schedules, dynamic values, and check for vaccine immunity reset.

Source code in CLT_BaseModel/flu_core/flu_components.py
def prepare_daily_state(self) -> None:
    """
    Override parent method to add vaccine immunity reset check.
    At beginning of each day, update schedules, dynamic values,
    and check for vaccine immunity reset.
    """
    # Call parent implementation first to update schedules and dynamic vals
    super().prepare_daily_state()

    # Check and potentially reset vaccine-induced immunity
    if hasattr(self.epi_metrics, 'MV'):
        self.epi_metrics.MV.check_and_apply_reset(
            self.current_real_date,
            self.params
        )

reset_simulation() -> None

Extends the base reset_simulation to recompute MV.init_val and vax_induced_*_risk_reduce_initial from the currently loaded vaccine schedule (and current base params) before resetting.

This ensures that if the daily_vaccines schedule has been replaced (e.g. via replace_schedule), or the underlying vax_induced_*_risk_reduce/vax_induced_immune_wane params have been overridden (e.g. by ScenarioRunner), the model resets to values consistent with the current schedule/params, rather than the values computed at construction time from the original schedule/params.

The MV.init_val recomputation uses VaxInducedImmunity.adjust_initial_value() with MV.original_init_val as the base β€” the unmodified value from the state JSON β€” so adjustments do not compound across calls.

Source code in CLT_BaseModel/flu_core/flu_components.py
def reset_simulation(self) -> None:
    """
    Extends the base `reset_simulation` to recompute `MV.init_val` and
    `vax_induced_*_risk_reduce_initial` from the currently loaded
    vaccine schedule (and current base params) before resetting.

    This ensures that if the `daily_vaccines` schedule has been replaced
    (e.g. via `replace_schedule`), or the underlying
    `vax_induced_*_risk_reduce`/`vax_induced_immune_wane` params have
    been overridden (e.g. by `ScenarioRunner`), the model resets to
    values consistent with the current schedule/params, rather than
    the values computed at construction time from the original
    schedule/params.

    The `MV.init_val` recomputation uses
    `VaxInducedImmunity.adjust_initial_value()` with
    `MV.original_init_val` as the base β€” the unmodified value
    from the state JSON β€” so adjustments do not compound across calls.
    """

    MV = self.epi_metrics["MV"]
    new_init_val = MV.adjust_initial_value(
        MV.original_init_val,
        self.start_real_date,
        self.params,
        self.schedules,
        self.simulation_settings.timesteps_per_day,
    )
    # Use the init_val setter so current_val is also updated immediately,
    # before super()'s reset loop overwrites it again (harmlessly).
    MV.init_val = new_init_val

    self.update_vax_induced_risk_reduce_initial()

    super().reset_simulation()

run_input_checks() -> None

Check the following
  • if total vaccinations exceed 100% over a year we issue a warning
  • school and work calendar variables must be between 0 and 1
  • absolute humidity values, contact matrix entries, daily vaccination must be non-negative
  • total contact matrix must be greater than the sum of the school and work contact matrices
  • all rate values must be strictly positive
  • initial compartmental values must be non-negative
Source code in CLT_BaseModel/flu_core/flu_components.py
def run_input_checks(self) -> None:
    """
    Check the following:
        - if total vaccinations exceed 100% over a year we issue a warning
        - school and work calendar variables must be between 0 and 1
        - absolute humidity values, contact matrix entries, daily vaccination
          must be non-negative
        - total contact matrix must be greater than the sum of the school
          and work contact matrices
        - all rate values must be strictly positive
        - initial compartmental values must be non-negative
    """

    self.check_humidity_input()
    self.check_vaccination_input()
    self.check_calendar_variables_input()
    self.check_contact_matrix_input()
    self.check_rate_input()
    self.check_initial_compartment_input()

update_vax_induced_risk_reduce_initial() -> None

Recomputes vax_induced_inf_risk_reduce_initial, vax_induced_hosp_risk_reduce_initial, and vax_induced_death_risk_reduce_initial from the current daily_vaccines schedule and updates self.params in place.

This must be re-run (not just computed once at construction) whenever the daily_vaccines schedule or any of the underlying base parameters (vax_induced_*_risk_reduce, vax_induced_immune_wane) change after construction -- e.g. via replace_schedule or a ScenarioRunner parameter override -- otherwise these derived values would silently keep reflecting the schedule/params from construction time. See reset_simulation, which calls this for the same reason VaxInducedImmunity's initial value is recomputed there.

Source code in CLT_BaseModel/flu_core/flu_components.py
def update_vax_induced_risk_reduce_initial(self) -> None:
    """
    Recomputes `vax_induced_inf_risk_reduce_initial`,
    `vax_induced_hosp_risk_reduce_initial`, and
    `vax_induced_death_risk_reduce_initial` from the current
    `daily_vaccines` schedule and updates `self.params` in place.

    This must be re-run (not just computed once at construction)
    whenever the `daily_vaccines` schedule or any of the underlying
    base parameters (`vax_induced_*_risk_reduce`,
    `vax_induced_immune_wane`) change after construction -- e.g. via
    `replace_schedule` or a `ScenarioRunner` parameter override --
    otherwise these derived values would silently keep reflecting
    the schedule/params from construction time. See
    `reset_simulation`, which calls this for the same reason
    `VaxInducedImmunity`'s initial value is recomputed there.
    """

    inf_initial, hosp_initial, death_initial = compute_vax_induced_risk_reduce_initial(
        self.params, self.schedules, self.start_real_date)
    self.params = clt.updated_dataclass(self.params, {
        "vax_induced_inf_risk_reduce_initial": inf_initial,
        "vax_induced_hosp_risk_reduce_initial": hosp_initial,
        "vax_induced_death_risk_reduce_initial": death_initial,
    })

FluSubpopModelError

Bases: SubpopModelError

Custom exceptions for flu subpopulation simulation model errors.

Source code in CLT_BaseModel/flu_core/flu_components.py
class FluSubpopModelError(clt.SubpopModelError):
    """Custom exceptions for flu subpopulation simulation model errors."""
    pass

FluSubpopParams dataclass

Bases: SubpopParams

Data container for pre-specified and fixed epidemiological parameters in FluSubpopModel.

Each field of datatype np.ndarray must be A x R, where A is the number of age groups and R is the number of risk groups. Note: this means all arrays should be 2D. See FluSubpopState docstring for important formatting note on 2D arrays.

Note: the user does not have to specify total_pop_age_risk -- this is automatically computed when a FluSubpopModel is instantiated. This is to ensure that the total population (summed across all compartments) actually equals total_pop_age_risk -- and the user doesn't change one without updating the other.

Attributes:

Name Type Description
num_age_groups positive int

number of age groups.

num_risk_groups positive int

number of risk groups.

start_real_date date

real-world date that corresponds to start sof simulation.

beta_baseline positive float

transmission rate.

total_pop_age_risk np.ndarray of positive ints

total number in population, summed across all age-risk groups.

humidity_impact positive float

coefficient that determines how much absolute humidity affects beta_baseline.

inf_induced_saturation np.ndarray of positive floats

constant(s) modeling saturation of antibody production of infected individuals.

inf_induced_immune_wane positive float

rate at which infection-induced immunity against infection wanes.

vax_induced_saturation np.ndarray of positive floats

constant(s) modeling saturation of antibody production of vaccinated individuals.

vax_induced_immune_wane positive float

rate at which vaccine-induced immunity against infection wanes.

inf_induced_inf_risk_reduce positive float

reduction in risk of getting infected after getting infected

inf_induced_hosp_risk_reduce positive float

reduction in risk of hospitalization after getting infected

inf_induced_death_risk_reduce positive float

reduction in risk of death after getting infected

vax_induced_inf_risk_reduce positive float

reduction in risk of getting infected after getting vaccinated

vax_induced_hosp_risk_reduce positive float

reduction in risk of hospitalization after getting vaccinated

vax_induced_death_risk_reduce positive float

reduction in risk of death after getting vaccinated

vax_induced_inf_risk_reduce_initial np.ndarray of positive floats

"peak" (zero-waning) vaccine-induced reduction in risk of getting infected -- computed once from vax_induced_inf_risk_reduce, vax_induced_immune_wane, and the daily_vaccines schedule so that, averaged over the vaccination season and accounting for waning, it reproduces vax_induced_inf_risk_reduce. See FluSubpopModel.compute_vax_induced_risk_reduce_initial.

vax_induced_hosp_risk_reduce_initial np.ndarray of positive floats

analogous "peak" (zero-waning) value for vax_induced_hosp_risk_reduce.

vax_induced_death_risk_reduce_initial np.ndarray of positive floats

analogous "peak" (zero-waning) value for vax_induced_death_risk_reduce.

adjust_VE_for_seasonal_waning bool

if True (default), vax_induced_inf_risk_reduce_initial, vax_induced_hosp_risk_reduce_initial, and vax_induced_death_risk_reduce_initial are computed from the corresponding vax_induced_*_risk_reduce season-average value, vax_induced_immune_wane, and the daily_vaccines schedule -- see compute_vax_induced_risk_reduce_initial. If False, this adjustment is skipped and the _initial values are simply set equal to the corresponding vax_induced_*_risk_reduce values (i.e. waning is not accounted for when determining the peak vaccine efficacy applied to newly-vaccinated individuals).

vax_protection_delay_days Optional[int]

(positive int): number of days after vaccination until vaccine protection is effective.

vax_immunity_reset_date_mm_dd Optional[str]

(str or None): date (in "mm_dd" format) each year when vaccine immunity resets, and date from which to start calculating contribution of vaccines to vaccine-induced immunity.

R_to_S_rate positive float

rate at which people in R move to S.

E_to_I_rate positive float

rate at which people in E move to I (both IP and IA, infected pre-symptomatic and infected asymptomatic)

IP_to_IS_rate positive float

rate a which people in IP (infected pre-symptomatic) move to IS (infected symptomatic)

ISR_to_R_rate positive float

rate at which people in IS (infected symptomatic) move to R.

IA_to_R_rate positive float

rate at which people in IA (infected asymptomatic) move to R

ISH_to_H_rate positive float

rate at which people in IS (infected symptomatic) move to H.

HR_to_R_rate positive float

rate at which people in H move to R.

HD_to_D_rate positive float

rate at which people in H move to D.

E_to_IA_prop np.ndarray of positive floats in [0,1]

proportion exposed who are asymptomatic based on age-risk groups.

IP_to_ISH_prop np.ndarray of positive floats in [0,1]

proportion infected who are hospitalized based on age-risk groups.

ISH_to_HD_prop np.ndarray of positive floats in [0,1]

proportion hospitalized who die based on age-risk groups.

IP_relative_inf positive float

relative infectiousness of pre-symptomatic to symptomatic people (IP to IS compartment).

IA_relative_inf positive float

relative infectiousness of asymptomatic to symptomatic people (IA to IS compartment).

relative_suscept np.ndarray of positive floats in [0,1]

relative susceptibility to infection by age group

total_contact_matrix np.ndarray of positive floats

A x A contact matrix (where A is the number of age groups), where element i,j is the average contacts from age group j that an individual in age group i has

school_contact_matrix np.ndarray of positive floats

A x A contact matrix (where A is the number of age groups), where element i,j is the average contacts from age group j that an individual in age group i has at school -- this matrix plus the work_contact_matrix must be less than the total_contact_matrix, element-wise

work_contact_matrix np.ndarray of positive floats

A x A contact matrix (where A is the number of age groups), where element i,j is the average contacts from age group j that an individual in age group i has at work -- this matrix plus the work_contact_matrix must be less than the total_contact_matrix, element-wise

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass(frozen=True)
class FluSubpopParams(clt.SubpopParams):
    """
    Data container for pre-specified and fixed epidemiological
    parameters in `FluSubpopModel`.

    Each field of datatype np.ndarray must be A x R,
    where A is the number of age groups and R is the number of
    risk groups. Note: this means all arrays should be 2D.
    See FluSubpopState docstring for important formatting note
    on 2D arrays.

    Note: the user does not have to specify `total_pop_age_risk` --
    this is automatically computed when a `FluSubpopModel` is
    instantiated. This is to ensure that the total population
    (summed across all compartments) actually equals `total_pop_age_risk` --
    and the user doesn't change one without updating the other.

    Attributes:
        num_age_groups (positive int):
            number of age groups.
        num_risk_groups (positive int):
            number of risk groups.
        start_real_date (datetime.date):
            real-world date that corresponds to start sof
            simulation.
        beta_baseline (positive float): transmission rate.
        total_pop_age_risk (np.ndarray of positive ints):
            total number in population, summed across all
            age-risk groups.
        humidity_impact (positive float):
            coefficient that determines how much absolute
            humidity affects beta_baseline.
        inf_induced_saturation (np.ndarray of positive floats):
            constant(s) modeling saturation of antibody
            production of infected individuals.
        inf_induced_immune_wane (positive float):
            rate at which infection-induced immunity
            against infection wanes.
        vax_induced_saturation (np.ndarray of positive floats):
            constant(s) modeling saturation of antibody
            production of vaccinated individuals.
        vax_induced_immune_wane (positive float):
            rate at which vaccine-induced immunity
            against infection wanes.
        inf_induced_inf_risk_reduce (positive float):
            reduction in risk of getting infected
            after getting infected
        inf_induced_hosp_risk_reduce (positive float):
            reduction in risk of hospitalization
            after getting infected
        inf_induced_death_risk_reduce (positive float):
            reduction in risk of death
            after getting infected
        vax_induced_inf_risk_reduce (positive float):
            reduction in risk of getting infected
            after getting vaccinated
        vax_induced_hosp_risk_reduce (positive float):
            reduction in risk of hospitalization
            after getting vaccinated
        vax_induced_death_risk_reduce (positive float):
            reduction in risk of death
            after getting vaccinated
        vax_induced_inf_risk_reduce_initial (np.ndarray of positive floats):
            "peak" (zero-waning) vaccine-induced reduction in risk
            of getting infected -- computed once from
            vax_induced_inf_risk_reduce, vax_induced_immune_wane,
            and the daily_vaccines schedule so that, averaged over
            the vaccination season and accounting for waning, it
            reproduces vax_induced_inf_risk_reduce. See
            FluSubpopModel.compute_vax_induced_risk_reduce_initial.
        vax_induced_hosp_risk_reduce_initial (np.ndarray of positive floats):
            analogous "peak" (zero-waning) value for
            vax_induced_hosp_risk_reduce.
        vax_induced_death_risk_reduce_initial (np.ndarray of positive floats):
            analogous "peak" (zero-waning) value for
            vax_induced_death_risk_reduce.
        adjust_VE_for_seasonal_waning (bool):
            if True (default), `vax_induced_inf_risk_reduce_initial`,
            `vax_induced_hosp_risk_reduce_initial`, and
            `vax_induced_death_risk_reduce_initial` are computed from
            the corresponding `vax_induced_*_risk_reduce` season-average
            value, `vax_induced_immune_wane`, and the `daily_vaccines`
            schedule -- see
            `compute_vax_induced_risk_reduce_initial`. If False, this
            adjustment is skipped and the `_initial` values are simply
            set equal to the corresponding `vax_induced_*_risk_reduce`
            values (i.e. waning is not accounted for when determining
            the peak vaccine efficacy applied to newly-vaccinated
            individuals).
        vax_protection_delay_days: (positive int):
            number of days after vaccination until vaccine
            protection is effective.
        vax_immunity_reset_date_mm_dd: (str or None):
            date (in "mm_dd" format) each year when vaccine
            immunity resets, and date from which to start
            calculating contribution of vaccines to 
            vaccine-induced immunity.
        R_to_S_rate (positive float):
            rate at which people in R move to S.
        E_to_I_rate (positive float):
            rate at which people in E move to I (both
            IP and IA, infected pre-symptomatic and infected
            asymptomatic)
        IP_to_IS_rate (positive float):
            rate a which people in IP (infected pre-symptomatic)
            move to IS (infected symptomatic)
        ISR_to_R_rate (positive float):
            rate at which people in IS (infected symptomatic)
            move to R.
        IA_to_R_rate (positive float):
            rate at which people in IA (infected asymptomatic)
            move to R
        ISH_to_H_rate (positive float):
            rate at which people in IS (infected symptomatic)
            move to H.
        HR_to_R_rate (positive float):
            rate at which people in H move to R.
        HD_to_D_rate (positive float):
            rate at which people in H move to D.
        E_to_IA_prop (np.ndarray of positive floats in [0,1]):
            proportion exposed who are asymptomatic based on
            age-risk groups.
        IP_to_ISH_prop (np.ndarray of positive floats in [0,1]):
            proportion infected who are hospitalized
            based on age-risk groups.
        ISH_to_HD_prop (np.ndarray of positive floats in [0,1]):
            proportion hospitalized who die based on
            age-risk groups.
        IP_relative_inf (positive float):
            relative infectiousness of pre-symptomatic to symptomatic
            people (IP to IS compartment).
        IA_relative_inf (positive float):
            relative infectiousness of asymptomatic to symptomatic
            people (IA to IS compartment).
        relative_suscept (np.ndarray of positive floats in [0,1]):
            relative susceptibility to infection by age group
        total_contact_matrix (np.ndarray of positive floats):
            A x A contact matrix (where A is the number
            of age groups), where element i,j is the average
            contacts from age group j that an individual in
            age group i has
        school_contact_matrix (np.ndarray of positive floats):
            A x A contact matrix (where A is the number
            of age groups), where element i,j is the average
            contacts from age group j that an individual in
            age group i has at school -- this matrix plus the
            work_contact_matrix must be less than the
            total_contact_matrix, element-wise
        work_contact_matrix (np.ndarray of positive floats):
            A x A contact matrix (where A is the number
            of age groups), where element i,j is the average
            contacts from age group j that an individual in
            age group i has at work -- this matrix plus the
            work_contact_matrix must be less than the
            total_contact_matrix, element-wise
    """

    num_age_groups: Optional[int] = None
    num_risk_groups: Optional[int] = None
    start_real_date: Optional[datetime.date] = None
    beta_baseline: Optional[float] = None
    total_pop_age_risk: Optional[np.ndarray] = None
    humidity_impact: Optional[float] = None

    inf_induced_saturation: Optional[float] = None
    inf_induced_immune_wane: Optional[float] = None
    vax_induced_saturation: Optional[float] = None
    vax_induced_immune_wane: Optional[float] = None
    inf_induced_inf_risk_reduce: Optional[float] = None
    inf_induced_hosp_risk_reduce: Optional[float] = None
    inf_induced_death_risk_reduce: Optional[float] = None
    vax_induced_inf_risk_reduce: Optional[float] = None
    vax_induced_hosp_risk_reduce: Optional[float] = None
    vax_induced_death_risk_reduce: Optional[float] = None
    vax_induced_inf_risk_reduce_initial: Optional[np.ndarray] = None
    vax_induced_hosp_risk_reduce_initial: Optional[np.ndarray] = None
    vax_induced_death_risk_reduce_initial: Optional[np.ndarray] = None
    adjust_VE_for_seasonal_waning: Optional[bool] = True
    vax_protection_delay_days: Optional[int] = 0
    vax_immunity_reset_date_mm_dd: Optional[str] = None

    R_to_S_rate: Optional[float] = None
    E_to_I_rate: Optional[float] = None
    IP_to_IS_rate: Optional[float] = None
    ISR_to_R_rate: Optional[float] = None
    IA_to_R_rate: Optional[float] = None
    ISH_to_H_rate: Optional[float] = None
    HR_to_R_rate: Optional[float] = None
    HD_to_D_rate: Optional[float] = None

    E_to_IA_prop: Optional[np.ndarray] = None
    IP_to_ISH_prop: Optional[torch.Tensor] = None
    ISH_to_HD_prop: Optional[torch.Tensor] = None

    IP_relative_inf: Optional[float] = None
    IA_relative_inf: Optional[float] = None

    relative_suscept: Optional[np.ndarray] = None

    total_contact_matrix: Optional[np.ndarray] = None
    school_contact_matrix: Optional[np.ndarray] = None
    work_contact_matrix: Optional[np.ndarray] = None

FluSubpopSchedules dataclass

Data container for dataframes used to specify schedules for each FluSubpopModel instance.

THE FORMAT FOR EACH DATAFRAME IS VERY IMPORTANT -- please read and implement carefully.

Attributes:

Name Type Description
absolute_humidity DataFrame

must have columns "date" and "absolute_humidity" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date objects -- "value" entries correspond to absolute humidity on those days

flu_contact_matrix DataFrame

must have columns "date", "is_school_day", and "is_work_day" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date object and "is_school_day" and "is_work_day" entries are are floats between 0 and 1 indicating if that date is a school day or work day

daily_vaccines DataFrame

must have "date" and "daily_vaccines" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date objects -- "value" entries correspond to historical proportion vaccinated on those days

mobility_modifier DataFrame

must have columns "mobility_modifier" and either "date" or "day_of_week" -- "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date objects -- "day_of_week" entries are strings with values from Monday to Sunday (case doesn't matter). "mobility_modifier" entries are JSON-encoded A x R arrays representing the proportion of time spent away from home by age-risk group on those days

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluSubpopSchedules:
    """
    Data container for dataframes used to specify schedules
    for each `FluSubpopModel` instance.

    THE FORMAT FOR EACH DATAFRAME IS VERY IMPORTANT -- please
    read and implement carefully.

    Attributes:
        absolute_humidity (pd.DataFrame):
            must have columns "date" and "absolute_humidity" --
            "date" entries must correspond to consecutive calendar days
            and must either be strings with `"YYYY-MM-DD"` format or
            `datetime.date` objects -- "value" entries correspond to
            absolute humidity on those days
        flu_contact_matrix (pd.DataFrame):
            must have columns "date", "is_school_day", and "is_work_day"
            -- "date" entries must correspond to consecutive calendar
            days and must either be strings with `"YYYY-MM-DD"` format
            or `datetime.date` object and "is_school_day" and
            "is_work_day" entries are are floats between 0 and 1 
            indicating if that date is a school day or work day
        daily_vaccines (pd.DataFrame):
            must have "date" and "daily_vaccines" -- "date" entries must
            correspond to consecutive calendar days and must either
            be strings with `"YYYY-MM-DD"` format or `datetime.date`
            objects -- "value" entries correspond to historical
            proportion vaccinated on those days
        mobility_modifier (pd.DataFrame):
            must have columns "mobility_modifier" and either "date"
            or "day_of_week" -- "date" entries must
            correspond to consecutive calendar days and must either
            be strings with `"YYYY-MM-DD"` format or `datetime.date`
            objects -- "day_of_week" entries are strings with values
            from Monday to Sunday (case doesn't matter).
            "mobility_modifier" entries are JSON-encoded A x R
            arrays representing the proportion of time spent away from
            home by age-risk group on those days
    """

    absolute_humidity: Optional[pd.DataFrame] = None
    flu_contact_matrix: Optional[pd.DataFrame] = None
    daily_vaccines: Optional[pd.DataFrame] = None
    mobility_modifier: Optional[pd.DataFrame] = None

FluSubpopState dataclass

Bases: SubpopState

Data container for pre-specified and fixed set of Compartment initial values and EpiMetric initial values for FluSubpopModel.

Each field below should be A x R np.ndarray, where A is the number of age groups and R is the number of risk groups. Note: this means all arrays should be 2D. Even if there is 1 age group and 1 risk group (no group stratification), each array should be 1x1, which is two-dimensional. For example, np.array([[100]]) is correct -- np.array([100]) is wrong.

Attributes:

Name Type Description
S np.ndarray of nonnegative integers

susceptible compartment for age-risk groups -- (holds current_val of Compartment "S").

E np.ndarray of nonnegative integers

exposed compartment for age-risk groups -- (holds current_val of Compartment "E").

IP np.ndarray of nonnegative integers

infected pre-symptomatic compartment for age-risk groups (holds current_val of Compartment "IP").

ISR np.ndarray of nonnegative integers

infected symptomatic (that will recover) compartment for age-risk groups (holds current_val of Compartment "ISR").

ISH np.ndarray of nonnegative integers

infected symptomatic compartment (that will be hospitalized) for age-risk groups (holds current_val of Compartment "ISH").

IA np.ndarray of nonnegative integers

infected asymptomatic compartment for age-risk groups (holds current_val of Compartment "IA").

HR np.ndarray of nonnegative integers

hospital compartment (that will recover) for age-risk groups (holds current_val of Compartment "HR").

HD np.ndarray of nonnegative integers

hospital compartment (that will die) for age-risk groups (holds current_val of Compartment "HD").

R np.ndarray of nonnegative integers

recovered compartment for age-risk groups (holds current_val of Compartment "R").

D np.ndarray of nonnegative integers

dead compartment for age-risk groups (holds current_val of Compartment "D").

M np.ndarray of nonnegative floats

infection-induced population-level immunity for age-risk groups (holds current_val of EpiMetric "M").

MV np.ndarray of nonnegative floats

vaccine-induced population-level immunity for age-risk groups (holds current_val of EpiMetric "MV").

absolute_humidity positive float

grams of water vapor per cubic meter g/m^3, used as seasonality parameter that influences transmission rate beta_baseline.

flu_contact_matrix np.ndarray of positive floats

A x A array, where A is the number of age groups -- element (a, a') corresponds to the number of contacts that a person in age group a has with people in age-risk group a'.

beta_reduce float in [0, 1]

starting value of DynamicVal "beta_reduce" on starting day of simulation -- this DynamicVal emulates a simple staged-alert policy

daily_vaccines np.ndarray of positive ints

holds current value of DailyVaccines instance, corresponding proportion of individuals in each age and risk group who received influenza vaccine on that day (generally derived from historical data)

mobility_modifier np.ndarray of positive floats

holds current value of MobilityModifier schedule, representing the proportion of time spent away from home by age group (A x R array)

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluSubpopState(clt.SubpopState):
    """
    Data container for pre-specified and fixed set of
    Compartment initial values and EpiMetric initial values
    for `FluSubpopModel`.

    Each field below should be A x R np.ndarray, where
    A is the number of age groups and R is the number of risk groups.
    Note: this means all arrays should be 2D. Even if there is
    1 age group and 1 risk group (no group stratification),
    each array should be 1x1, which is two-dimensional.
    For example, np.array([[100]]) is correct --
    np.array([100]) is wrong.

    Attributes:
        S (np.ndarray of nonnegative integers):
            susceptible compartment for age-risk groups --
            (holds current_val of Compartment "S").
        E (np.ndarray of nonnegative integers):
            exposed compartment for age-risk groups --
            (holds current_val of Compartment "E").
        IP (np.ndarray of nonnegative integers):
            infected pre-symptomatic compartment for age-risk groups
            (holds current_val of Compartment "IP").
        ISR (np.ndarray of nonnegative integers):
            infected symptomatic (that will recover) compartment
            for age-risk groups
            (holds current_val of Compartment "ISR").
        ISH (np.ndarray of nonnegative integers):
            infected symptomatic compartment (that will be hospitalized)
            for age-risk groups
            (holds current_val of Compartment "ISH").
        IA (np.ndarray of nonnegative integers):
            infected asymptomatic compartment for age-risk groups
            (holds current_val of Compartment "IA").
        HR (np.ndarray of nonnegative integers):
            hospital compartment (that will recover)
            for age-risk groups
            (holds current_val of Compartment "HR").
        HD (np.ndarray of nonnegative integers):
            hospital compartment (that will die)
            for age-risk groups
            (holds current_val of Compartment "HD").
        R (np.ndarray of nonnegative integers):
            recovered compartment for age-risk groups
            (holds current_val of Compartment "R").
        D (np.ndarray of nonnegative integers):
            dead compartment for age-risk groups
            (holds current_val of Compartment "D").
        M (np.ndarray of nonnegative floats):
            infection-induced population-level immunity
            for age-risk groups (holds current_val
            of EpiMetric "M").
        MV (np.ndarray of nonnegative floats):
            vaccine-induced population-level immunity
            for age-risk groups (holds current_val
            of EpiMetric "MV").
        absolute_humidity (positive float):
            grams of water vapor per cubic meter g/m^3,
            used as seasonality parameter that influences
            transmission rate beta_baseline.
        flu_contact_matrix (np.ndarray of positive floats):
            A x A array, where A is the number of age
            groups -- element (a, a') corresponds to the number
            of contacts that a person in age group a
            has with people in age-risk group a'.
        beta_reduce (float in [0,1]):
            starting value of DynamicVal "beta_reduce" on
            starting day of simulation -- this DynamicVal
            emulates a simple staged-alert policy
        daily_vaccines (np.ndarray of positive ints):
            holds current value of DailyVaccines instance,
            corresponding proportion of individuals in each
            age and risk group who received influenza
            vaccine on that day (generally derived from
            historical data)
        mobility_modifier (np.ndarray of positive floats):
            holds current value of MobilityModifier schedule,
            representing the proportion of time spent away
            from home by age group (A x R array)
    """

    S: Optional[np.ndarray] = None
    E: Optional[np.ndarray] = None
    IP: Optional[np.ndarray] = None
    ISR: Optional[np.ndarray] = None
    ISH: Optional[np.ndarray] = None
    IA: Optional[np.ndarray] = None
    HR: Optional[np.ndarray] = None
    HD: Optional[np.ndarray] = None
    R: Optional[np.ndarray] = None
    D: Optional[np.ndarray] = None

    M: Optional[np.ndarray] = None
    MV: Optional[np.ndarray] = None

    absolute_humidity: Optional[float] = None
    flu_contact_matrix: Optional[np.ndarray] = None
    beta_reduce: Optional[float] = 0.0

    daily_vaccines: Optional[np.ndarray] = None
    mobility_modifier: Optional[np.ndarray] = None

FluTravelParamsTensors dataclass

Data container for tensors for FluMetapopModel -- used to store arrays that contain data across all subpopulations (collected from parameters on each location/subpopulation model, as well as from the metapopulation's associated FluMixingParams instance). Note that not all fields in FluSubpopParams are included -- we only include parameters needed for the travel model computation, for efficiency.

Attributes:

Name Type Description
num_locations (Tensor, 0 - dimensional)

number of locations (subpopulations) in the metapopulation model and therefore the travel model.

travel_proportions Tensor

L x L array, where L is the number of locations or subpopulations, where element i,j corresponds to the proportion of the population in location i who travels to location j (on average).

See FluSubpopParams docstring for other attributes.

Fields are analogous -- but (most) are size (L, A, R) for location-age-risk or size 0 tensors. Exceptions are travel_proportions, which is size (L, L), and any of the contact matrices, which are size (L, A, A).

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluTravelParamsTensors:
    """
    Data container for tensors for `FluMetapopModel` -- used to store arrays
    that contain data across all subpopulations (collected from parameters
    on each location/subpopulation model, as well as from the
    metapopulation's associated `FluMixingParams` instance).
    Note that not all fields in `FluSubpopParams` are included
    -- we only include parameters needed for the travel model
    computation, for efficiency.

    Attributes:
        num_locations (torch.Tensor, 0-dimensional):
            number of locations (subpopulations) in the
            metapopulation model and therefore the travel
            model.
        travel_proportions (torch.Tensor):
            L x L array, where L is the number of locations
            or subpopulations, where element i,j corresponds
            to the proportion of the population in location i
            who travels to location j (on average).

    See `FluSubpopParams` docstring for other attributes.

    Fields are analogous -- but (most) are size (L, A, R) for
    location-age-risk or size 0 tensors. Exceptions are
    `travel_proportions`, which is size (L, L),
    and any of the contact matrices, which are size (L, A, A).
    """

    num_locations: Optional[torch.tensor] = None
    num_age_groups: Optional[torch.tensor] = None
    num_risk_groups: Optional[torch.tensor] = None

    travel_proportions: torch.Tensor = None

    IP_relative_inf: torch.Tensor = None
    IA_relative_inf: torch.Tensor = None

    relative_suscept: torch.Tensor = None

    total_contact_matrix: Optional[torch.Tensor] = None
    school_contact_matrix: Optional[torch.Tensor] = None
    work_contact_matrix: Optional[torch.Tensor] = None

    def standardize_shapes(self) -> None:
        """
        If field is size (L, A, R) for location-age-risk or size 0 tensors,
            or is not a special variable listed below, then apply dimension
            expansion so that fields are size (L, A, R) tensors for tensor multiplication.

        Exceptions are `travel_proportions`, which is size (L, L),
        and any of the contact matrices, which are size (L, A, A).

        Not all dimension combinations are considered not all make sense --
        we assume that we only have risk IF we have age, for example.
        """

        L = int(self.num_locations.item())
        A = int(self.num_age_groups.item())
        R = int(self.num_risk_groups.item())

        error_str = " Each SubpopParams field must have size (L, A, R) " \
                    "(for location-age-risk groups) or size 0 -- please check files " \
                    "and inputs, then try again."

        for name, value in vars(self).items():

            # Ignore the field that corresponds to a dictionary
            if name == "init_vals":
                continue

            elif name == "travel_proportions":
                if value.size() != torch.Size([L, L]):
                    raise Exception(str(name) + error_str)

            # `total_contact_matrix`, `school_contact_matrix`, `work_contact_matrix`
            elif "contact_matrix" in name:
                if value.size() == torch.Size([L, A, A]):
                    continue
                elif value.size() != torch.Size([A, A]):
                    raise Exception(str(name) + error_str)
                else:
                    setattr(self, name, value.view(1, A, A).expand(L, A, A))

            # string parameters
            elif isinstance(value, str) or isinstance(value, datetime.date):
                continue

            # If scalar or already L x A x R, do not need to adjust
            #   dimensions
            elif value.size() == torch.Size([]):
                continue

            elif value.size() == torch.Size([L, A, R]):
                continue

            elif value.size() == torch.Size([L]):
                setattr(self, name, value.view(L, 1, 1).expand(L, A, R))

            elif value.size() == torch.Size([A, R]):
                setattr(self, name, value.view(1, A, R).expand(L, A, R))

            else:
                value_size = str(value.size())
                raise Exception(str(name) + ' with size ' + value_size + error_str)

standardize_shapes() -> None

If field is size (L, A, R) for location-age-risk or size 0 tensors, or is not a special variable listed below, then apply dimension expansion so that fields are size (L, A, R) tensors for tensor multiplication.

Exceptions are travel_proportions, which is size (L, L), and any of the contact matrices, which are size (L, A, A).

Not all dimension combinations are considered not all make sense -- we assume that we only have risk IF we have age, for example.

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
def standardize_shapes(self) -> None:
    """
    If field is size (L, A, R) for location-age-risk or size 0 tensors,
        or is not a special variable listed below, then apply dimension
        expansion so that fields are size (L, A, R) tensors for tensor multiplication.

    Exceptions are `travel_proportions`, which is size (L, L),
    and any of the contact matrices, which are size (L, A, A).

    Not all dimension combinations are considered not all make sense --
    we assume that we only have risk IF we have age, for example.
    """

    L = int(self.num_locations.item())
    A = int(self.num_age_groups.item())
    R = int(self.num_risk_groups.item())

    error_str = " Each SubpopParams field must have size (L, A, R) " \
                "(for location-age-risk groups) or size 0 -- please check files " \
                "and inputs, then try again."

    for name, value in vars(self).items():

        # Ignore the field that corresponds to a dictionary
        if name == "init_vals":
            continue

        elif name == "travel_proportions":
            if value.size() != torch.Size([L, L]):
                raise Exception(str(name) + error_str)

        # `total_contact_matrix`, `school_contact_matrix`, `work_contact_matrix`
        elif "contact_matrix" in name:
            if value.size() == torch.Size([L, A, A]):
                continue
            elif value.size() != torch.Size([A, A]):
                raise Exception(str(name) + error_str)
            else:
                setattr(self, name, value.view(1, A, A).expand(L, A, A))

        # string parameters
        elif isinstance(value, str) or isinstance(value, datetime.date):
            continue

        # If scalar or already L x A x R, do not need to adjust
        #   dimensions
        elif value.size() == torch.Size([]):
            continue

        elif value.size() == torch.Size([L, A, R]):
            continue

        elif value.size() == torch.Size([L]):
            setattr(self, name, value.view(L, 1, 1).expand(L, A, R))

        elif value.size() == torch.Size([A, R]):
            setattr(self, name, value.view(1, A, R).expand(L, A, R))

        else:
            value_size = str(value.size())
            raise Exception(str(name) + ' with size ' + value_size + error_str)

FluTravelStateTensors dataclass

Data container for tensors for FluMetapopModel -- used to store arrays that contain data across all subpopulations (collected from each location/subpopulation model). Note that not all fields in FluSubpopState are included -- we only include compartments needed for the travel model computation, for efficiency.

Attributes:

Name Type Description
IP torch.Tensor of nonnegative integers

presymptomatic infected compartment for location-age-risk groups -- the lth element holds current_val of Compartment "IP" on the lth location / subpopulation on the associated MetapopModel.

ISR torch.Tensor of nonnegative integers

symptomatic infected compartment (that will recover) for location-age-risk groups -- the lth element holds current_val of Compartment "ISR" on the lth location / subpopulation on the associated MetapopModel.

ISH torch.Tensor of nonnegative integers

symptomatic infected compartment (that will be hospitalized) for location-age-risk groups -- the lth element holds current_val of Compartment "ISH" on the lth location / subpopulation on the associated MetapopModel.

IA torch.Tensor of nonnegative integers

asymptomatic infected compartment for location-age-risk groups -- the lth element holds current_val of Compartment "IA" on the lth location / subpopulation on the associated MetapopModel.

HR torch.Tensor of nonnegative integers

hospital compartment (that will recover) for location-age-risk groups -- the lth element holds current_val of Compartment "HR" on the lth location / subpopulation on the associated MetapopModel.

HD torch.Tensor of nonnegative integers

hospital compartment (that will die) for location-age-risk groups -- the lth element holds current_val of Compartment "HD" on the lth location / subpopulation on the associated MetapopModel.

flu_contact_matrix torch.Tensor of nonnegative integers

contact matrix for location-age-risk groups -- the lth element holds current_val of FluContactMatrix Schedule for subpopulation l -- this value is a combination of the total contact matrix, the work contact matrix, and the school contact matrix (and the value is adjusted depending on whether the date is a work or school day)

init_vals dict

dictionary of torch.Tensor instances, where keys correspond to "IP", "ISR", "ISH", "IA", "HR", and "HD", and values correspond to their initial values for location-age-risk groups.

mobility_modifier torch.Tensor of positive floats

mobility modifier for location-age-risk groups -- the lth element holds current_val of MobilityModifier Schedule for subpopulation l -- represents the proportion of time spent away from home by age group

Source code in CLT_BaseModel/flu_core/flu_data_structures.py
@dataclass
class FluTravelStateTensors:
    """
    Data container for tensors for `FluMetapopModel` -- used to store arrays
    that contain data across all subpopulations (collected from each
    location/subpopulation model). Note that not all fields in
    `FluSubpopState` are included -- we only include compartments
    needed for the travel model computation, for efficiency.

    Attributes:
        IP (torch.Tensor of nonnegative integers):
            presymptomatic infected compartment for location-age-risk
            groups -- the lth element holds current_val of
            Compartment "IP" on the lth location / subpopulation
            on the associated `MetapopModel`.
        ISR (torch.Tensor of nonnegative integers):
            symptomatic infected compartment (that will recover)
            for location-age-risk groups -- the lth element holds
            current_val of Compartment "ISR" on the lth 
            location / subpopulation on the associated `MetapopModel`.
        ISH (torch.Tensor of nonnegative integers):
            symptomatic infected compartment (that will be hospitalized)
            for location-age-risk groups -- the lth element holds
            current_val of Compartment "ISH" on the lth 
            location / subpopulation on the associated `MetapopModel`.
        IA (torch.Tensor of nonnegative integers):
            asymptomatic infected compartment for location-age-risk
            groups -- the lth element holds current_val of
            Compartment "IA" on the lth location / subpopulation
            on the associated `MetapopModel`.
        HR (torch.Tensor of nonnegative integers):
            hospital compartment (that will recover) for location-age-risk
            groups -- the lth element holds current_val of
            Compartment "HR" on the lth location / subpopulation
            on the associated `MetapopModel`.
        HD (torch.Tensor of nonnegative integers):
            hospital compartment (that will die) for location-age-risk
            groups -- the lth element holds current_val of
            Compartment "HD" on the lth location / subpopulation
            on the associated `MetapopModel`.
        flu_contact_matrix (torch.Tensor of nonnegative integers):
            contact matrix for location-age-risk groups -- the
            lth element holds current_val of `FluContactMatrix`
            `Schedule` for subpopulation l -- this value is a
            combination of the total contact matrix, the
            work contact matrix, and the school contact matrix
            (and the value is adjusted depending on whether
            the date is a work or school day)
        init_vals (dict):
            dictionary of torch.Tensor instances, where keys
            correspond to "IP", "ISR", "ISH", "IA", "HR", and "HD", and
            values correspond to their initial values for
            location-age-risk groups.
        mobility_modifier (torch.Tensor of positive floats):
            mobility modifier for location-age-risk groups -- the
            lth element holds current_val of `MobilityModifier`
            `Schedule` for subpopulation l -- represents the proportion
            of time spent away from home by age group
    """

    IP: torch.Tensor = None
    ISR: torch.Tensor = None
    ISH: torch.Tensor = None
    IA: torch.Tensor = None
    HR: torch.Tensor = None
    HD: torch.Tensor = None

    flu_contact_matrix: torch.Tensor = None
    mobility_modifier: torch.Tensor = None

    init_vals: dict = field(default_factory=dict)

    # Note: `init_vals: dict = {}` does NOT work --
    #   gives "mutable default" argument

    def save_current_vals_as_init_vals(self):

        for field in fields(self):
            if field.name == "init_vals":
                continue
            self.init_vals[field.name] = getattr(self, field.name).clone()

    def reset_to_init_vals(self):

        for name, val in self.init_vals.items():
            setattr(self, name, val.clone())

HospDeadToDead

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "HD" to "D" compartment. The functional form is the same across subpopulations.

Source code in CLT_BaseModel/flu_core/flu_components.py
class HospDeadToDead(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "HD" to "D" compartment. The functional form is the same across
    subpopulations.
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.HD_to_D_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.HD_to_D_rate)

HospRecoverToRecovered

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "HR" to "R" compartment. The functional form is the same across subpopulations.

Source code in CLT_BaseModel/flu_core/flu_components.py
class HospRecoverToRecovered(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "HR" to "R" compartment. The functional form is the same across
    subpopulations.
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.HR_to_R_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.HR_to_R_rate)

InfInducedImmunity

Bases: EpiMetric

EpiMetric-derived class for infection-induced population-level immunity.

Population-level immunity increases as people move from "R" to "S" -- this is a design choice intended to avoid "double-counting." People in "R" cannot be infected at all. People who move from "R" to "S" are susceptible again, but these recently-recovered people should have partial immunity. To handle this phenomenon, this epi metric increases as people move from "R" to "S."

Parameters:

Name Type Description Default
R_to_S RecoveredToSusceptible

RecoveredToSusceptible TransitionVariable in the SubpopModel -- it is an attribute because the population-level immunity increases as people move from "R" to "S".

required

See parent class docstring for other attributes.

Source code in CLT_BaseModel/flu_core/flu_components.py
class InfInducedImmunity(clt.EpiMetric):
    """
    EpiMetric-derived class for infection-induced
    population-level immunity.

    Population-level immunity increases as people move
    from "R" to "S" -- this is a design choice intended
    to avoid "double-counting." People in "R" cannot be
    infected at all. People who move from "R" to "S"
    are susceptible again, but these recently-recovered people
    should have partial immunity. To handle this phenomenon,
    this epi metric increases as people move from "R" to "S."

    Params:
        R_to_S (RecoveredToSusceptible):
            RecoveredToSusceptible TransitionVariable
            in the SubpopModel -- it is an attribute
            because the population-level immunity
            increases as people move from "R" to "S".

    See parent class docstring for other attributes.
    """

    def __init__(self, init_val, R_to_S):
        super().__init__(init_val)
        self.R_to_S = R_to_S

    def get_change_in_current_val(self,
                                  state: FluSubpopState,
                                  params: FluSubpopParams,
                                  num_timesteps: int) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        # Note: the current values of transition variables already include
        #   discretization (division by the number of timesteps) -- therefore,
        #   we do not divide the first part of this equation by the number of
        #   timesteps -- see `TransitionVariable` class's methods for getting
        #   various realizations for more information

        return (self.R_to_S.current_val / params.total_pop_age_risk) * \
               (1 - params.inf_induced_saturation * state.M) - \
               params.inf_induced_immune_wane * state.M / num_timesteps

get_change_in_current_val(state: FluSubpopState, params: FluSubpopParams, num_timesteps: int) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_change_in_current_val(self,
                              state: FluSubpopState,
                              params: FluSubpopParams,
                              num_timesteps: int) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    # Note: the current values of transition variables already include
    #   discretization (division by the number of timesteps) -- therefore,
    #   we do not divide the first part of this equation by the number of
    #   timesteps -- see `TransitionVariable` class's methods for getting
    #   various realizations for more information

    return (self.R_to_S.current_val / params.total_pop_age_risk) * \
           (1 - params.inf_induced_saturation * state.M) - \
           params.inf_induced_immune_wane * state.M / num_timesteps

MobilityModifier

Bases: Schedule

Schedule for time-varying mobility modifier values.

Attributes:

Name Type Description
timeseries_df DataFrame

There are 2 possible input formats: i) a standard schedule that must have columns "date" and "mobility_modifier" where "date" entries must correspond to consecutive calendar days and must either be strings with "YYYY-MM-DD" format or datetime.date objects ii) a day of week schedule that must have columns "day_of_week" and "mobility_modifier" where "day_of_week" entries are strings with values from Monday to Sunday (case doesn't matter). The code will automatically detect which format is being used by looking at the column name. In both cases, "mobility_modifier" entries are JSON-encoded A x R arrays representing the proportion of time spent away from home by age-risk group on those days. Identical to FluSubpopSchedules field of same name.

Source code in CLT_BaseModel/flu_core/flu_components.py
class MobilityModifier(clt.Schedule):
    """
    Schedule for time-varying mobility modifier values.

    Attributes:
        timeseries_df (pd.DataFrame):
            There are 2 possible input formats:
            i) a standard schedule that must have columns "date" and
            "mobility_modifier" where "date" entries must correspond
            to consecutive calendar days and must either be strings with
            `"YYYY-MM-DD"` format or `datetime.date` objects 
            ii) a day of week schedule that must have columns "day_of_week"
            and "mobility_modifier" where "day_of_week" entries are
            strings with values from Monday to Sunday (case doesn't matter).
            The code will automatically detect which format is being used
            by looking at the column name.
            In both cases, "mobility_modifier" entries are
            JSON-encoded A x R arrays representing the proportion of
            time spent away from home by age-risk group on those days.
            Identical to `FluSubpopSchedules` field of same name.
    """

    def __init__(self,
                 init_val: Optional[np.ndarray | float] = None,
                 timeseries_df: pd.DataFrame = None):
        """
        Args:
            init_val (Optional[np.ndarray | float]):
                starting value(s) at the beginning of the simulation
            timeseries_df (Optional[pd.DataFrame] = None):
                must have columns ("date" or "day_of_week") 
                and "mobility_modifier" --
                see class docstring for format details.
        """

        super().__init__(init_val)

        self.timeseries_df = timeseries_df

    def update_current_val(self, params, current_date: datetime.date) -> None:
        if self.is_day_of_week_schedule:
            current_day_of_week = current_date.strftime('%A').lower()
            self.current_val = self.timeseries_df.loc[current_day_of_week, "mobility_modifier"]
        else:
            self.current_val = self.timeseries_df.loc[current_date, "mobility_modifier"]

    def postprocess_data_input(self) -> None:
        """
            Converts mobility_modifier column from
            a string representation of a list of lists
            (each day) of format AxR into np.ndarray.
            Check whether day_of_week schedule is being used.
            Make days of week lower case if being used.
            Pre-indexes the DataFrame by date or day_of_week for O(1) lookups.
        """

        if 'day_of_week' in self.timeseries_df.columns:
            self.is_day_of_week_schedule = True

        self.timeseries_df['mobility_modifier'] = \
            self.timeseries_df['mobility_modifier'].apply(json.loads)
        self.timeseries_df.loc[:, 'mobility_modifier'] = \
            self.timeseries_df['mobility_modifier'].apply(
                lambda x: np.asarray(x)
                )

        if self.is_day_of_week_schedule:
            self.timeseries_df['day_of_week'] = \
                self.timeseries_df['day_of_week'].str.lower()
            self.timeseries_df = self.timeseries_df.set_index('day_of_week')
        else:
            self.timeseries_df = self.timeseries_df.set_index('date')

__init__(init_val: Optional[np.ndarray | float] = None, timeseries_df: pd.DataFrame = None)

Parameters:

Name Type Description Default
init_val Optional[ndarray | float]

starting value(s) at the beginning of the simulation

None
timeseries_df Optional[pd.DataFrame] = None

must have columns ("date" or "day_of_week") and "mobility_modifier" -- see class docstring for format details.

None
Source code in CLT_BaseModel/flu_core/flu_components.py
def __init__(self,
             init_val: Optional[np.ndarray | float] = None,
             timeseries_df: pd.DataFrame = None):
    """
    Args:
        init_val (Optional[np.ndarray | float]):
            starting value(s) at the beginning of the simulation
        timeseries_df (Optional[pd.DataFrame] = None):
            must have columns ("date" or "day_of_week") 
            and "mobility_modifier" --
            see class docstring for format details.
    """

    super().__init__(init_val)

    self.timeseries_df = timeseries_df

postprocess_data_input() -> None

Converts mobility_modifier column from a string representation of a list of lists (each day) of format AxR into np.ndarray. Check whether day_of_week schedule is being used. Make days of week lower case if being used. Pre-indexes the DataFrame by date or day_of_week for O(1) lookups.

Source code in CLT_BaseModel/flu_core/flu_components.py
def postprocess_data_input(self) -> None:
    """
        Converts mobility_modifier column from
        a string representation of a list of lists
        (each day) of format AxR into np.ndarray.
        Check whether day_of_week schedule is being used.
        Make days of week lower case if being used.
        Pre-indexes the DataFrame by date or day_of_week for O(1) lookups.
    """

    if 'day_of_week' in self.timeseries_df.columns:
        self.is_day_of_week_schedule = True

    self.timeseries_df['mobility_modifier'] = \
        self.timeseries_df['mobility_modifier'].apply(json.loads)
    self.timeseries_df.loc[:, 'mobility_modifier'] = \
        self.timeseries_df['mobility_modifier'].apply(
            lambda x: np.asarray(x)
            )

    if self.is_day_of_week_schedule:
        self.timeseries_df['day_of_week'] = \
            self.timeseries_df['day_of_week'].str.lower()
        self.timeseries_df = self.timeseries_df.set_index('day_of_week')
    else:
        self.timeseries_df = self.timeseries_df.set_index('date')

PresympToSympHospital

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "IP" to "ISH" compartment. The functional form is the same across subpopulations.

Each PresympToSympHospital instance forms a TransitionVariableGroup with a corresponding PresympToSympRecover instance (these two transition variables are jointly distributed).

Source code in CLT_BaseModel/flu_core/flu_components.py
class PresympToSympHospital(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "IP" to "ISH" compartment. The functional form is the same across
    subpopulations.

    Each PresympToSympHospital instance forms a TransitionVariableGroup with
    a corresponding PresympToSympRecover instance (these two
    transition variables are jointly distributed).
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """
        inf_induced_hosp_risk_reduce = params.inf_induced_hosp_risk_reduce
        inf_induced_proportional_risk_reduce = inf_induced_hosp_risk_reduce / (1 - inf_induced_hosp_risk_reduce)

        immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

        vax_immunity_factor = 1 - state.MV * params.vax_induced_hosp_risk_reduce_initial

        prob_hosp = (params.IP_to_ISH_prop / immunity_force) * vax_immunity_factor

        return np.asarray(prob_hosp * params.IP_to_IS_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """
    inf_induced_hosp_risk_reduce = params.inf_induced_hosp_risk_reduce
    inf_induced_proportional_risk_reduce = inf_induced_hosp_risk_reduce / (1 - inf_induced_hosp_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_hosp_risk_reduce_initial

    prob_hosp = (params.IP_to_ISH_prop / immunity_force) * vax_immunity_factor

    return np.asarray(prob_hosp * params.IP_to_IS_rate)

PresympToSympRecover

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "IP" to "ISR" compartment. The functional form is the same across subpopulations.

Each PresympToSympRecover instance forms a TransitionVariableGroup with a corresponding PresympToSympHospital instance (these two transition variables are jointly distributed).

Source code in CLT_BaseModel/flu_core/flu_components.py
class PresympToSympRecover(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "IP" to "ISR" compartment. The functional form is the same across
    subpopulations.

    Each PresympToSympRecover instance forms a TransitionVariableGroup with
    a corresponding PresympToSympHospital instance (these two
    transition variables are jointly distributed).
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """
        inf_induced_hosp_risk_reduce = params.inf_induced_hosp_risk_reduce
        inf_induced_proportional_risk_reduce = inf_induced_hosp_risk_reduce / (1 - inf_induced_hosp_risk_reduce)

        immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

        vax_immunity_factor = 1 - state.MV * params.vax_induced_hosp_risk_reduce_initial

        prob_hosp = (params.IP_to_ISH_prop / immunity_force) * vax_immunity_factor

        return np.asarray((1 - prob_hosp) * params.IP_to_IS_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """
    inf_induced_hosp_risk_reduce = params.inf_induced_hosp_risk_reduce
    inf_induced_proportional_risk_reduce = inf_induced_hosp_risk_reduce / (1 - inf_induced_hosp_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_hosp_risk_reduce_initial

    prob_hosp = (params.IP_to_ISH_prop / immunity_force) * vax_immunity_factor

    return np.asarray((1 - prob_hosp) * params.IP_to_IS_rate)

RecoveredToSusceptible

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "R" to "S" compartment. The functional form is the same across subpopulations.

Source code in CLT_BaseModel/flu_core/flu_components.py
class RecoveredToSusceptible(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "R" to "S" compartment. The functional form is the same across
    subpopulations.
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """
        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.R_to_S_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """
    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.R_to_S_rate)

SusceptibleToExposed

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "S" to "E" compartment. The functional form is the same across subpopulations.

The rate depends on the corresponding subpopulation's contact matrix, transmission rate beta, number infected (symptomatic, asymptomatic, and pre-symptomatic), and population-level immunity against infection, among other parameters.

This is the most complicated transition variable in the flu model. If using metapopulation model (travel model), then the rate depends on the total_mixing_exposure attribute, which is a function of other subpopulations' states and parameters, and travel between subpopulations.

If there is no metapopulation model, the rate is much simpler.

Attributes:

Name Type Description
total_mixing_exposure np.ndarray of positive floats

weighted infectious count (exposure) from movement within home location, travel to other locations, and visitors from other locations

See parent class docstring for other attributes.

Source code in CLT_BaseModel/flu_core/flu_components.py
class SusceptibleToExposed(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "S" to "E" compartment. The functional form is the same across
    subpopulations.

    The rate depends on the corresponding subpopulation's
    contact matrix, transmission rate beta, number
    infected (symptomatic, asymptomatic, and pre-symptomatic),
    and population-level immunity against infection,
    among other parameters.

    This is the most complicated transition variable in the
    flu model. If using metapopulation model (travel model), then
    the rate depends on the `total_mixing_exposure` attribute,
    which is a function of other subpopulations' states and
    parameters, and travel between subpopulations.

    If there is no metapopulation model, the rate
    is much simpler.

    Attributes:
        total_mixing_exposure (np.ndarray of positive floats):
            weighted infectious count (exposure) from movement
            within home location, travel to other locations,
            and visitors from other locations

    See parent class docstring for other attributes.
    """

    def __init__(self,
                 origin: clt.Compartment,
                 destination: clt.Compartment,
                 transition_type: clt.TransitionTypes,
                 is_jointly_distributed: str = False):

        super().__init__(origin,
                         destination,
                         transition_type,
                         is_jointly_distributed)

        self.total_mixing_exposure = None

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        # If `total_mixing_exposure` has not been updated,
        #   then there is no travel model -- so, simulate
        #   this subpopulation entirely independently and
        #   use the simplified transition rate that does not
        #   depend on travel dynamics

        beta_adjusted = compute_beta_adjusted(state, params)

        inf_induced_inf_risk_reduce = params.inf_induced_inf_risk_reduce
        inf_induced_proportional_risk_reduce = inf_induced_inf_risk_reduce / (1 - inf_induced_inf_risk_reduce)

        immune_force = 1 + inf_induced_proportional_risk_reduce * state.M

        # Vaccine-induced protection against infection is modeled as a
        #   multiplicative reduction (rather than folded into the additive
        #   `immune_force` denominator above) -- see
        #   `compute_vax_induced_risk_reduce_initial` for how
        #   `vax_induced_inf_risk_reduce_initial` is derived.
        vax_immunity_factor = 1 - state.MV * params.vax_induced_inf_risk_reduce_initial

        if self.total_mixing_exposure is not None:

            # Note here `self.total_mixing_exposure` includes
            #   `suscept_by_age` -- see `compute_total_mixing_exposure_prop`
            #   in `flu_travel_functions`

            # Need to convert tensor into array because combining np.ndarrays and
            #   tensors doesn't work, and everything else is an array
            # Note: `self.total_mixing_exposure` (a Tensor) must be the left
            #   operand of the first multiplication -- np.ndarray * Tensor
            #   raises TypeError, but Tensor * np.ndarray works fine.
            return np.asarray(
                (beta_adjusted * self.total_mixing_exposure / immune_force) * vax_immunity_factor)

        else:
            wtd_presymp_asymp_by_age = compute_wtd_presymp_asymp_by_age(state, params)

            # Super confusing syntax... but this is the pain of having A x R,
            #   but having the contact matrix (contact patterns) be for
            #   ONLY age groups
            wtd_infectious_prop = np.divide(np.sum(sum([state.ISR, state.ISH]), axis=1, keepdims=True) + wtd_presymp_asymp_by_age, # Remy TODO check the sum works
                                            compute_pop_by_age(params))

            raw_total_exposure = np.matmul(state.flu_contact_matrix, wtd_infectious_prop)

            # The total rate is only age-dependent -- it's the same rate across age groups
            return params.relative_suscept * (beta_adjusted * vax_immunity_factor * raw_total_exposure / immune_force)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    # If `total_mixing_exposure` has not been updated,
    #   then there is no travel model -- so, simulate
    #   this subpopulation entirely independently and
    #   use the simplified transition rate that does not
    #   depend on travel dynamics

    beta_adjusted = compute_beta_adjusted(state, params)

    inf_induced_inf_risk_reduce = params.inf_induced_inf_risk_reduce
    inf_induced_proportional_risk_reduce = inf_induced_inf_risk_reduce / (1 - inf_induced_inf_risk_reduce)

    immune_force = 1 + inf_induced_proportional_risk_reduce * state.M

    # Vaccine-induced protection against infection is modeled as a
    #   multiplicative reduction (rather than folded into the additive
    #   `immune_force` denominator above) -- see
    #   `compute_vax_induced_risk_reduce_initial` for how
    #   `vax_induced_inf_risk_reduce_initial` is derived.
    vax_immunity_factor = 1 - state.MV * params.vax_induced_inf_risk_reduce_initial

    if self.total_mixing_exposure is not None:

        # Note here `self.total_mixing_exposure` includes
        #   `suscept_by_age` -- see `compute_total_mixing_exposure_prop`
        #   in `flu_travel_functions`

        # Need to convert tensor into array because combining np.ndarrays and
        #   tensors doesn't work, and everything else is an array
        # Note: `self.total_mixing_exposure` (a Tensor) must be the left
        #   operand of the first multiplication -- np.ndarray * Tensor
        #   raises TypeError, but Tensor * np.ndarray works fine.
        return np.asarray(
            (beta_adjusted * self.total_mixing_exposure / immune_force) * vax_immunity_factor)

    else:
        wtd_presymp_asymp_by_age = compute_wtd_presymp_asymp_by_age(state, params)

        # Super confusing syntax... but this is the pain of having A x R,
        #   but having the contact matrix (contact patterns) be for
        #   ONLY age groups
        wtd_infectious_prop = np.divide(np.sum(sum([state.ISR, state.ISH]), axis=1, keepdims=True) + wtd_presymp_asymp_by_age, # Remy TODO check the sum works
                                        compute_pop_by_age(params))

        raw_total_exposure = np.matmul(state.flu_contact_matrix, wtd_infectious_prop)

        # The total rate is only age-dependent -- it's the same rate across age groups
        return params.relative_suscept * (beta_adjusted * vax_immunity_factor * raw_total_exposure / immune_force)

SympHospitalToHospDead

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "ISH" to "HD" compartment. The functional form is the same across subpopulations.

Each SympHospitalToHospDead instance forms a TransitionVariableGroup with a corresponding SympHospitalToHospRecover instance (these two transition variables are jointly distributed).

The rate of SympHospitalToHospDead decreases as population-level immunity against hospitalization increases.

Source code in CLT_BaseModel/flu_core/flu_components.py
class SympHospitalToHospDead(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "ISH" to "HD" compartment. The functional form is the same across
    subpopulations.

    Each SympHospitalToHospDead instance forms a TransitionVariableGroup with
    a corresponding SympHospitalToHospRecover instance (these two
    transition variables are jointly distributed).

    The rate of SympHospitalToHospDead decreases as population-level immunity
    against hospitalization increases.
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        inf_induced_death_risk_reduce = params.inf_induced_death_risk_reduce

        inf_induced_proportional_risk_reduce = \
            inf_induced_death_risk_reduce / (1 - inf_induced_death_risk_reduce)

        immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

        vax_immunity_factor = 1 - state.MV * params.vax_induced_death_risk_reduce_initial

        prob_death = (params.ISH_to_HD_prop / immunity_force) * vax_immunity_factor

        return np.asarray(prob_death * params.ISH_to_H_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    inf_induced_death_risk_reduce = params.inf_induced_death_risk_reduce

    inf_induced_proportional_risk_reduce = \
        inf_induced_death_risk_reduce / (1 - inf_induced_death_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_death_risk_reduce_initial

    prob_death = (params.ISH_to_HD_prop / immunity_force) * vax_immunity_factor

    return np.asarray(prob_death * params.ISH_to_H_rate)

SympHospitalToHospRecover

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "ISH" to "HR" compartment. The functional form is the same across subpopulations.

Each SympHospitalToHospRecover instance forms a TransitionVariableGroup with a corresponding SympHospitalToHospDead instance (these two transition variables are jointly distributed).

Source code in CLT_BaseModel/flu_core/flu_components.py
class SympHospitalToHospRecover(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "ISH" to "HR" compartment. The functional form is the same across
    subpopulations.

    Each SympHospitalToHospRecover instance forms a TransitionVariableGroup with
    a corresponding SympHospitalToHospDead instance (these two
    transition variables are jointly distributed).
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        inf_induced_death_risk_reduce = params.inf_induced_death_risk_reduce

        inf_induced_proportional_risk_reduce = \
            inf_induced_death_risk_reduce / (1 - inf_induced_death_risk_reduce)

        immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

        vax_immunity_factor = 1 - state.MV * params.vax_induced_death_risk_reduce_initial

        prob_death = (params.ISH_to_HD_prop / immunity_force) * vax_immunity_factor

        return np.asarray((1 - prob_death) * params.ISH_to_H_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    inf_induced_death_risk_reduce = params.inf_induced_death_risk_reduce

    inf_induced_proportional_risk_reduce = \
        inf_induced_death_risk_reduce / (1 - inf_induced_death_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_death_risk_reduce_initial

    prob_death = (params.ISH_to_HD_prop / immunity_force) * vax_immunity_factor

    return np.asarray((1 - prob_death) * params.ISH_to_H_rate)

SympRecoverToRecovered

Bases: TransitionVariable

TransitionVariable-derived class for movement from the "ISR" to "R" compartment. The functional form is the same across subpopulations.

Source code in CLT_BaseModel/flu_core/flu_components.py
class SympRecoverToRecovered(clt.TransitionVariable):
    """
    TransitionVariable-derived class for movement from the
    "ISR" to "R" compartment. The functional form is the same across
    subpopulations.
    """

    def get_current_rate(self,
                         state: FluSubpopState,
                         params: FluSubpopParams) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        return np.full((params.num_age_groups, params.num_risk_groups),
                       params.ISR_to_R_rate)

get_current_rate(state: FluSubpopState, params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_current_rate(self,
                     state: FluSubpopState,
                     params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    return np.full((params.num_age_groups, params.num_risk_groups),
                   params.ISR_to_R_rate)

VaxInducedImmunity

Bases: EpiMetric

EpiMetric-derived class for vaccine-induced population-level immunity.

Source code in CLT_BaseModel/flu_core/flu_components.py
class VaxInducedImmunity(clt.EpiMetric):
    """
    EpiMetric-derived class for vaccine-induced
    population-level immunity.
    """

    def __init__(self, 
                 init_val,
                 current_real_date: datetime.date,
                 params: FluSubpopParams,
                 schedules: clt.Schedule,
                 timesteps_per_day: int):

        adjusted_init_val = self.adjust_initial_value(
            init_val, current_real_date, params, schedules, timesteps_per_day)
        super().__init__(adjusted_init_val)


    def adjust_initial_value(self,
                             init_val: np.ndarray,
                             current_real_date: datetime.date,
                             params: FluSubpopParams,
                             schedules: clt.Schedule,
                             timesteps_per_day: int):
        """
        Adjusts initial value of vaccine-induced immunity based on
        vaccination schedule if the parameter vax_immunity_reset_date_mm_dd
        is not None.
        Vaccines administered before the reset date are not counted, but
        vaccines administered after the reset date (and before the simulation
        start date, accounting for protection delay) are counted with waning.
        """

        self.original_init_val = copy.deepcopy(init_val)
        self.adjusted_init_val = copy.deepcopy(init_val)

        if params.vax_immunity_reset_date_mm_dd is not None:
            # Print warning to mention that initial value is being adjusted
            msg = 'Vaccine immunity reset date is set as ' +\
                f'{params.vax_immunity_reset_date_mm_dd.replace("_", "/")}. \n' +\
                'Initial vaccine-induced immunity value is being adjusted ' +\
                'by resetting immunity to 0 at that date, and by taking into ' +\
                'account vaccines administered after this date, and before simulation start date.'
            warnings.warn(msg)

            # Parse reset date (format: "MM_DD")
            month, day = params.vax_immunity_reset_date_mm_dd.split('_')
            current_year = current_real_date.year

            # Find most recent occurrence of reset date before start_date
            reset_date = datetime.date(current_year, int(month), int(day))
            if reset_date >= current_real_date:
                # If reset date is after start, use previous year
                reset_date = datetime.date(current_year - 1, int(month), int(day))

            # Filter vaccines between reset_date and start date,
            # accounting for protection delay
            vaccines_df = schedules['daily_vaccines'].timeseries_df.copy()

            mask = (vaccines_df.index >= (reset_date + datetime.timedelta(days=params.vax_protection_delay_days))) &\
                (vaccines_df.index < current_real_date)
            relevant_vaccines = vaccines_df[mask]

            # Initialize MV adjustment
            MV_adjustment = np.zeros_like(vaccines_df['daily_vaccines'].iloc[0])

            # Process each day's vaccines, applying waning
            for _, row in relevant_vaccines.iterrows():
                for i in range(timesteps_per_day):
                    MV_adjustment += row["daily_vaccines"] / timesteps_per_day - \
                        params.vax_induced_immune_wane * MV_adjustment / timesteps_per_day

            self.adjusted_init_val = self.adjusted_init_val + MV_adjustment

        return self.adjusted_init_val

    def get_change_in_current_val(self,
                                  state: FluSubpopState,
                                  params: FluSubpopParams,
                                  num_timesteps: int) -> np.ndarray:
        """
        Returns:
            np.ndarray of shape (A, R)
        """

        # Note: `state.daily_vaccines` (based on the value of the `DailyVaccines`
        #   `Schedule` is NOT divided by the number of timesteps -- so we need to
        #   do this division in the equation here.

        return state.daily_vaccines / (num_timesteps) - \
               params.vax_induced_immune_wane * state.MV / num_timesteps

    def check_and_apply_reset(self, 
                              current_date: datetime.date, 
                              params: FluSubpopParams):
        """
        Check if current date matches vax_immunity_reset_date_mm_dd.
        If so, reset MV to zero.

        Args:
            current_date: The current simulation date
            params: FluSubpopParams containing vax_immunity_reset_date_mm_dd
        """

        if params.vax_immunity_reset_date_mm_dd is not None:
            # Parse reset date (format: "MM_DD")
            month, day = params.vax_immunity_reset_date_mm_dd.split('_')

            # Check if current date matches the reset date (month and day)
            if current_date.month == int(month) and current_date.day == int(day):
                # Reset vaccine-induced immunity to zero
                self.current_val = np.zeros_like(self.current_val)
                print(f"VaxInducedImmunity reset to 0 on {current_date}")

adjust_initial_value(init_val: np.ndarray, current_real_date: datetime.date, params: FluSubpopParams, schedules: clt.Schedule, timesteps_per_day: int)

Adjusts initial value of vaccine-induced immunity based on vaccination schedule if the parameter vax_immunity_reset_date_mm_dd is not None. Vaccines administered before the reset date are not counted, but vaccines administered after the reset date (and before the simulation start date, accounting for protection delay) are counted with waning.

Source code in CLT_BaseModel/flu_core/flu_components.py
def adjust_initial_value(self,
                         init_val: np.ndarray,
                         current_real_date: datetime.date,
                         params: FluSubpopParams,
                         schedules: clt.Schedule,
                         timesteps_per_day: int):
    """
    Adjusts initial value of vaccine-induced immunity based on
    vaccination schedule if the parameter vax_immunity_reset_date_mm_dd
    is not None.
    Vaccines administered before the reset date are not counted, but
    vaccines administered after the reset date (and before the simulation
    start date, accounting for protection delay) are counted with waning.
    """

    self.original_init_val = copy.deepcopy(init_val)
    self.adjusted_init_val = copy.deepcopy(init_val)

    if params.vax_immunity_reset_date_mm_dd is not None:
        # Print warning to mention that initial value is being adjusted
        msg = 'Vaccine immunity reset date is set as ' +\
            f'{params.vax_immunity_reset_date_mm_dd.replace("_", "/")}. \n' +\
            'Initial vaccine-induced immunity value is being adjusted ' +\
            'by resetting immunity to 0 at that date, and by taking into ' +\
            'account vaccines administered after this date, and before simulation start date.'
        warnings.warn(msg)

        # Parse reset date (format: "MM_DD")
        month, day = params.vax_immunity_reset_date_mm_dd.split('_')
        current_year = current_real_date.year

        # Find most recent occurrence of reset date before start_date
        reset_date = datetime.date(current_year, int(month), int(day))
        if reset_date >= current_real_date:
            # If reset date is after start, use previous year
            reset_date = datetime.date(current_year - 1, int(month), int(day))

        # Filter vaccines between reset_date and start date,
        # accounting for protection delay
        vaccines_df = schedules['daily_vaccines'].timeseries_df.copy()

        mask = (vaccines_df.index >= (reset_date + datetime.timedelta(days=params.vax_protection_delay_days))) &\
            (vaccines_df.index < current_real_date)
        relevant_vaccines = vaccines_df[mask]

        # Initialize MV adjustment
        MV_adjustment = np.zeros_like(vaccines_df['daily_vaccines'].iloc[0])

        # Process each day's vaccines, applying waning
        for _, row in relevant_vaccines.iterrows():
            for i in range(timesteps_per_day):
                MV_adjustment += row["daily_vaccines"] / timesteps_per_day - \
                    params.vax_induced_immune_wane * MV_adjustment / timesteps_per_day

        self.adjusted_init_val = self.adjusted_init_val + MV_adjustment

    return self.adjusted_init_val

check_and_apply_reset(current_date: datetime.date, params: FluSubpopParams)

Check if current date matches vax_immunity_reset_date_mm_dd. If so, reset MV to zero.

Parameters:

Name Type Description Default
current_date date

The current simulation date

required
params FluSubpopParams

FluSubpopParams containing vax_immunity_reset_date_mm_dd

required
Source code in CLT_BaseModel/flu_core/flu_components.py
def check_and_apply_reset(self, 
                          current_date: datetime.date, 
                          params: FluSubpopParams):
    """
    Check if current date matches vax_immunity_reset_date_mm_dd.
    If so, reset MV to zero.

    Args:
        current_date: The current simulation date
        params: FluSubpopParams containing vax_immunity_reset_date_mm_dd
    """

    if params.vax_immunity_reset_date_mm_dd is not None:
        # Parse reset date (format: "MM_DD")
        month, day = params.vax_immunity_reset_date_mm_dd.split('_')

        # Check if current date matches the reset date (month and day)
        if current_date.month == int(month) and current_date.day == int(day):
            # Reset vaccine-induced immunity to zero
            self.current_val = np.zeros_like(self.current_val)
            print(f"VaxInducedImmunity reset to 0 on {current_date}")

get_change_in_current_val(state: FluSubpopState, params: FluSubpopParams, num_timesteps: int) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def get_change_in_current_val(self,
                              state: FluSubpopState,
                              params: FluSubpopParams,
                              num_timesteps: int) -> np.ndarray:
    """
    Returns:
        np.ndarray of shape (A, R)
    """

    # Note: `state.daily_vaccines` (based on the value of the `DailyVaccines`
    #   `Schedule` is NOT divided by the number of timesteps -- so we need to
    #   do this division in the equation here.

    return state.daily_vaccines / (num_timesteps) - \
           params.vax_induced_immune_wane * state.MV / num_timesteps

accept_reject_admits(metapop_model: FluMetapopModel, sampling_RNG: np.random.Generator, sampling_info: dict[str, dict[str, clt.UniformSamplingSpec]], total_daily_target_admits: list[np.ndarray], num_days: int = 50, target_accepted_reps: int = int(100.0), max_reps: int = int(1000.0), early_stop_percent: float = 0.5, target_rsquared: float = 0.75)

Accept-reject sampler for a metapopulation model.

This function repeatedly samples parameters from uniform distributions (as specified in spec) and simulates the model until the R-squared between simulated total admits and reference data exceeds target_rsquared. Accepted parameter sets and simulation states are saved as JSON files.

Parameters:

Name Type Description Default
metapop_model FluMetapopModel

The metapopulation model to simulate and sample parameters for.

required
sampling_RNG Generator

Random number generator used for uniform sampling.

required
sampling_info dict[str, dict[str, UniformSamplingSpec]]

See clt_toolkit / sampling / sample_uniform_metapop_params / sampling_info parameter for description.

required
total_daily_target_admits list[ndarray]

"Target" time series of total admits (across subpopulations) for computing R-squared -- we would like parameters and sample paths that give simulated admits close to total_daily_target_admits. Must have length equal to num_days.

required
num_days int, default=50

Total number of days to simulate for accepted parameter sets.

50
target_accepted_reps int, default=100

Target number of accepted parameter sets (replicates) to collect.

int(100.0)
max_reps int, default=1000

Maximum number of sampling attempts before stopping.

int(1000.0)
early_stop_percent float, default=0.5

Fraction of num_days to simulate initially for early R-squared check.

0.5
target_rsquared float, default=0.75

Minimum R-squared required between simulated and reference admits for acceptance.

0.75

Notes: - Early stopping is performed at num_days * early_stop_percent to quickly reject poor parameter samples. - Accepted samples (and the state of the simulation at day num_days) are saved to JSON files per subpopulation. Note that for efficiency, NOT ALL PARAMETERS ARE SAVED! Only the parameters that are randomly sampled (and thus are different between replications). - Running this function can be slow -- test this function with a small number of replications or simulation days to start.

Source code in CLT_BaseModel/flu_core/flu_accept_reject.py
def accept_reject_admits(metapop_model: FluMetapopModel,
                         sampling_RNG: np.random.Generator,
                         sampling_info: dict[str, dict[str, clt.UniformSamplingSpec]],
                         total_daily_target_admits: list[np.ndarray],
                         num_days: int = 50,
                         target_accepted_reps: int = int(1e2),
                         max_reps: int = int(1e3),
                         early_stop_percent: float = 0.5,
                         target_rsquared: float = 0.75):
    """
    Accept-reject sampler for a metapopulation model.

    This function repeatedly samples parameters from uniform distributions
    (as specified in `spec`) and simulates the model until the R-squared between
    simulated total admits and reference data exceeds `target_rsquared`.
    Accepted parameter sets and simulation states are saved as JSON files.

    Parameters:
        metapop_model (flu.FluMetapopModel):
            The metapopulation model to simulate and sample parameters for.
        sampling_RNG (np.random.Generator):
            Random number generator used for uniform sampling.
        sampling_info (dict[str, dict[str, clt.UniformSamplingSpec]]):
            See `clt_toolkit / sampling / sample_uniform_metapop_params / sampling_info`
            parameter for description.
        total_daily_target_admits (list[np.ndarray]):
            "Target" time series of total admits (across subpopulations)
            for computing R-squared -- we would like parameters and
            sample paths that give simulated admits close to
            `total_daily_target_admits`. Must have length equal to `num_days`.
        num_days (int, default=50):
            Total number of days to simulate for accepted parameter sets.
        target_accepted_reps (int, default=100):
            Target number of accepted parameter sets (replicates) to collect.
        max_reps (int, default=1000):
            Maximum number of sampling attempts before stopping.
        early_stop_percent (float, default=0.5):
            Fraction of `num_days` to simulate initially for early R-squared check.
        target_rsquared (float, default=0.75):
            Minimum R-squared required between simulated and reference admits for acceptance.

    Notes:
    - Early stopping is performed at `num_days * early_stop_percent` to
        quickly reject poor parameter samples.
    - Accepted samples (and the state of the simulation at day
        `num_days`) are saved to JSON files per subpopulation.
        Note that for efficiency, NOT ALL PARAMETERS ARE SAVED!
        Only the parameters that are randomly sampled (and thus are
        different between replications).
    - Running this function can be slow -- test this function with a small
        number of replications or simulation days to start.
    """

    if target_accepted_reps > max_reps:
        max_reps = 10 * target_accepted_reps

    num_days_early_stop = int(num_days * early_stop_percent)

    reps_counter = 0
    accepted_reps_counter = 0

    while reps_counter < max_reps and accepted_reps_counter < target_accepted_reps:

        reps_counter += 1

        metapop_model.reset_simulation()

        param_samples = clt.sample_uniform_metapop_params(metapop_model,
                                                          sampling_RNG,
                                                          sampling_info)

        # Save IS to H transition variable history
        # But do not save daily (compartment) history for efficiency
        for subpop_name, updates_dict in param_samples.items():
            metapop_model.modify_subpop_params(subpop_name, updates_dict)
            metapop_model.modify_simulation_settings({"transition_variables_to_save": ["ISH_to_HR", "ISH_to_HD"],
                                                      "save_daily_history": False})

        metapop_model.simulate_until_day(num_days_early_stop)
        total_simulated_admits = clt.aggregate_daily_tvar_history(metapop_model, ["ISH_to_HR", "ISH_to_HD"])
        current_rsquared = compute_rsquared(reference_timeseries=total_daily_target_admits[:num_days_early_stop],
                                            simulated_timeseries=total_simulated_admits)
        if current_rsquared < target_rsquared:
            continue

        else:
            metapop_model.simulate_until_day(num_days)
            total_simulated_admits = clt.aggregate_daily_tvar_history(metapop_model, ["ISH_to_HR", "ISH_to_HD"])
            current_rsquared = compute_rsquared(reference_timeseries=total_daily_target_admits,
                                                simulated_timeseries=total_simulated_admits)
            if current_rsquared < target_rsquared:
                continue
            else:
                accepted_reps_counter += 1

                for subpop_name, subpop in metapop_model.subpop_models.items():
                    with open("subpop_" + str(subpop_name) + "_rep_" + str(accepted_reps_counter) +
                              "_accepted_sample_params.json", "w") as f:
                        json.dump(clt.serialize_dataclass(param_samples[subpop_name]), f, indent=4)
                    with open("subpop_" + str(subpop_name) + "_rep_" + str(accepted_reps_counter) +
                              "_accepted_state.json", "w") as f:
                        json.dump(clt.serialize_dataclass(subpop.state), f, indent=4)

advance_timestep(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, precomputed: FluPrecomputedTensors, dt: float, save_calibration_targets: bool = False, save_tvar_history: bool = False, total_mixing_exposure: torch.Tensor = None) -> Tuple[FluFullMetapopStateTensors, dict, dict]

Advance the simulation one timestep, with length dt. Updates state corresponding to compartments and epidemiological metrics after computing transition variables and metric changes.

Note that in this torch "mean" deterministic implementation... - We compute rates in the same way as the get_binom_deterministic_no_round transition type in the OOP code -- see TransitionVariables class in clt_toolkit / base_components for more details. - We also implement a "mean" deterministic analog of the multinomial distribution to handle multiple outflows from the same compartment - We do not round the transition variables - We also use softplus, a smooth approximation to the ReLU function, to ensure that compartments are nonnegative (which is not guaranteed using the mean of a binomial/multinomial random variable rather than sampling from those distributions).

Returns:

Type Description
Tuple[FluFullMetapopStateTensors, dict, dict]

New FluFullMetapopStateTensors with updated state, dict of calibration targets corresponding to state values or transition variable values used for calibration, and dict of transition variable values to save this history. If save_calibration_targets is False, then the corresponding dict is empty, and similarly with save_tvar_history.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def advance_timestep(state: FluFullMetapopStateTensors,
                     params: FluFullMetapopParamsTensors,
                     precomputed: FluPrecomputedTensors,
                     dt: float,
                     save_calibration_targets: bool=False,
                     save_tvar_history: bool=False,
                     total_mixing_exposure: torch.Tensor = None) -> Tuple[FluFullMetapopStateTensors, dict, dict]:
    """
    Advance the simulation one timestep, with length `dt`.
    Updates state corresponding to compartments and
    epidemiological metrics after computing transition variables
    and metric changes.

    Note that in this torch "mean" deterministic implementation...
    - We compute rates in the same way as the
        `get_binom_deterministic_no_round`
        transition type in the OOP code -- see
        `TransitionVariables` class in
        `clt_toolkit / base_components` for more details.
    - We also implement a "mean" deterministic analog
        of the multinomial distribution to handle
        multiple outflows from the same compartment
    - We do not round the transition variables
    - We also use `softplus`, a smooth approximation to the
        ReLU function, to ensure that compartments are
        nonnegative (which is not guaranteed using
        the mean of a binomial/multinomial random variable
        rather than sampling from those distributions).

    Returns:
        (Tuple[FluFullMetapopStateTensors, dict, dict]):
            New `FluFullMetapopStateTensors` with updated state,
            `dict` of calibration targets corresponding to state
            values or transition variable values used for calibration,
            and `dict` of transition variable values to save this
            history. If `save_calibration_targets` is `False`,
            then the corresponding `dict` is empty, and similarly with
            `save_tvar_history`.
    """ 

    S_to_E = compute_S_to_E(state, params, precomputed, dt,
                            total_mixing_exposure=total_mixing_exposure)

    # Deterministic multinomial implementation to match
    #   object-oriented version
    E_to_IP_rate = compute_E_to_IP_rate(params)
    E_to_IA_rate = compute_E_to_IA_rate(params)
    E_outgoing_total_rate = E_to_IP_rate + E_to_IA_rate
    E_outgoing_total = state.E * \
        torch_approx_binom_probability_from_rate(E_outgoing_total_rate, dt)
    E_to_IA = E_outgoing_total * (E_to_IA_rate / E_outgoing_total_rate)              
    E_to_IP = E_outgoing_total * (E_to_IP_rate / E_outgoing_total_rate)

    IA_to_R = compute_IA_to_R(state, params, dt)

    # Deterministic multinomial implementation to match
    #   object-oriented version
    IP_to_ISR_rate = compute_IP_to_ISR_rate(state, params)
    IP_to_ISH_rate = compute_IP_to_ISH_rate(state, params)
    IP_outgoing_total_rate = IP_to_ISR_rate + IP_to_ISH_rate
    IP_outgoing_total = state.IP * \
        torch_approx_binom_probability_from_rate(IP_outgoing_total_rate, dt)
    IP_to_ISR = IP_outgoing_total * (IP_to_ISR_rate / IP_outgoing_total_rate)
    IP_to_ISH = IP_outgoing_total * (IP_to_ISH_rate / IP_outgoing_total_rate)

    ISR_to_R = compute_ISR_to_R(state, params, dt)

    # Deterministic multinomial implementation to match
    #   object-oriented version
    ISH_to_HR_rate = compute_ISH_to_HR_rate(state, params)
    ISH_to_HD_rate = compute_ISH_to_HD_rate(state, params)
    ISH_outgoing_total_rate = ISH_to_HR_rate + ISH_to_HD_rate
    ISH_outgoing_total = state.ISH * \
        torch_approx_binom_probability_from_rate(ISH_outgoing_total_rate, dt)
    ISH_to_HR = ISH_outgoing_total * (ISH_to_HR_rate / ISH_outgoing_total_rate)
    ISH_to_HD = ISH_outgoing_total * (ISH_to_HD_rate / ISH_outgoing_total_rate)

    # Deterministic multinomial implementation to match
    #   object-oriented version
    HR_to_R = compute_HR_to_R(state, params, dt)
    HD_to_D = compute_HD_to_D(state, params, dt)

    R_to_S = compute_R_to_S(state, params, dt)


    # Make sure compartments are nonnegative
    S_new = torch.nn.functional.softplus(state.S + R_to_S - S_to_E)
    E_new = torch.nn.functional.softplus(state.E + S_to_E - E_to_IP - E_to_IA)
    IP_new = torch.nn.functional.softplus(state.IP + E_to_IP - IP_to_ISR - IP_to_ISH)
    ISR_new = torch.nn.functional.softplus(state.ISR + IP_to_ISR - ISR_to_R)
    ISH_new = torch.nn.functional.softplus(state.ISH + IP_to_ISH - ISH_to_HR - ISH_to_HD)
    IA_new = torch.nn.functional.softplus(state.IA + E_to_IA - IA_to_R)
    HR_new = torch.nn.functional.softplus(state.HR + ISH_to_HR - HR_to_R)
    HD_new = torch.nn.functional.softplus(state.HD + ISH_to_HD - HD_to_D)
    R_new = torch.nn.functional.softplus(state.R + ISR_to_R + IA_to_R + HR_to_R - R_to_S)
    D_new = torch.nn.functional.softplus(state.D + HD_to_D)

    # Update immunity variables
    M_change = compute_M_change(state, params, precomputed, dt)
    MV_change = compute_MV_change(state, params, precomputed, dt)
    M_new = state.M + M_change
    MV_new = state.MV + MV_change

    state_new = FluFullMetapopStateTensors(S=S_new,
                                           E=E_new,
                                           IP=IP_new,
                                           ISR=ISR_new,
                                           ISH=ISH_new,
                                           IA=IA_new,
                                           HR=HR_new,
                                           HD=HD_new,
                                           R=R_new,
                                           D=D_new,
                                           M=M_new,
                                           MV=MV_new,
                                           absolute_humidity=state.absolute_humidity,
                                           daily_vaccines=state.daily_vaccines,
                                           flu_contact_matrix=state.flu_contact_matrix,
                                           mobility_modifier=state.mobility_modifier)

    calibration_targets = {}
    if save_calibration_targets:
        calibration_targets["ISH_to_H"] = ISH_to_HR + ISH_to_HD

    transition_variables = {}
    if save_tvar_history:
        transition_variables["S_to_E"] = S_to_E
        transition_variables["E_to_IP"] = E_to_IP
        transition_variables["E_to_IA"] = E_to_IA
        transition_variables["IA_to_R"] = IA_to_R
        transition_variables["IP_to_ISR"] = IP_to_ISR
        transition_variables["IP_to_ISH"] = IP_to_ISH
        transition_variables["ISR_to_R"] = ISR_to_R
        transition_variables["ISH_to_HR"] = ISH_to_HR
        transition_variables["ISH_to_HD"] = ISH_to_HD
        transition_variables["HR_to_R"] = HR_to_R
        transition_variables["HD_to_D"] = HD_to_D
        transition_variables["R_to_S"] = R_to_S
        transition_variables["M_change"] = M_change
        transition_variables["MV_change"] = MV_change

    return state_new, calibration_targets, transition_variables

attack_rate(metapop_model: MetapopModel, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None) -> float

Attack rate = cumulative infections / initial susceptible population.

Returns

float

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def attack_rate(
    metapop_model: MetapopModel,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
) -> float:
    """
    Attack rate = cumulative infections / initial susceptible population.

    Returns
    -------
    float
    """
    infections = daily_new_infections(
        metapop_model, subpop_name, age_group, risk_group
    ).sum()

    if subpop_name is not None:
        subpops = [metapop_model.subpop_models[subpop_name]]
    else:
        subpops = list(metapop_model.subpop_models.values())

    # Initial susceptible = S compartment at first recorded timestep
    init_S_arrays = [
        np.asarray(subpop.compartments["S"].history_vals_list[0])
        for subpop in subpops
    ]
    init_S = np.sum(np.stack(init_S_arrays, axis=0), axis=0)  # (A, R)

    if age_group is not None:
        init_S = init_S[age_group : age_group + 1, :]
    if risk_group is not None:
        init_S = init_S[:, risk_group : risk_group + 1]

    return float(infections / init_S.sum())

check_and_apply_MV_reset(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, day_counter: int)

Check if current date matches vax_immunity_reset_date_mm_dd where vaccine-induced immunity should be reset. If so, reset MV to zero.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def check_and_apply_MV_reset(state: FluFullMetapopStateTensors,
                             params: FluFullMetapopParamsTensors,
                             day_counter: int):
    """
    Check if current date matches vax_immunity_reset_date_mm_dd
    where vaccine-induced immunity should be reset.
    If so, reset MV to zero.
    """

    if params.vax_immunity_reset_date_mm_dd is not None:
        current_date = params.start_real_date + datetime.timedelta(days=day_counter)

        # Parse reset date (format: "MM_DD")
        month, day = params.vax_immunity_reset_date_mm_dd.split('_')

        # Check if current date matches the reset date (month and day)
        if current_date.month == int(month) and current_date.day == int(day):
            # Reset vaccine-induced immunity to zero
            state.MV = np.zeros_like(state.MV)
            print(f"VaxInducedImmunity MV reset to 0 on {current_date}")

compute_E_to_IA_rate(params: FluFullMetapopParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_E_to_IA_rate(params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    return params.E_to_I_rate * params.E_to_IA_prop

compute_E_to_IP_rate(params: FluFullMetapopParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_E_to_IP_rate(params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    return params.E_to_I_rate * (1 - params.E_to_IA_prop)

compute_HD_to_D(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_HD_to_D(state: FluFullMetapopStateTensors,
                    params: FluFullMetapopParamsTensors,
                    dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    rate = params.HD_to_D_rate

    HD_to_D = state.HD * torch_approx_binom_probability_from_rate(rate, dt)

    return HD_to_D

compute_HR_to_R(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_HR_to_R(state: FluFullMetapopStateTensors,
                    params: FluFullMetapopParamsTensors,
                    dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    rate = params.HR_to_R_rate

    HR_to_R = state.HR * torch_approx_binom_probability_from_rate(rate, dt)

    return HR_to_R

compute_IA_to_R(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_IA_to_R(state: FluFullMetapopStateTensors,
                    params: FluFullMetapopParamsTensors,
                    dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    rate = params.IA_to_R_rate

    IA_to_R = state.IA * torch_approx_binom_probability_from_rate(rate, dt)

    return IA_to_R

compute_IP_to_ISH_rate(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_IP_to_ISH_rate(state: FluFullMetapopStateTensors,
                           params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    inf_induced_hosp_risk_reduce = params.inf_induced_hosp_risk_reduce
    inf_induced_proportional_risk_reduce = inf_induced_hosp_risk_reduce / (1 - inf_induced_hosp_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_hosp_risk_reduce_initial

    prob_hosp = (params.IP_to_ISH_prop / immunity_force) * vax_immunity_factor

    rate = params.IP_to_IS_rate * prob_hosp

    return rate

compute_IP_to_ISR_rate(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_IP_to_ISR_rate(state: FluFullMetapopStateTensors,
                           params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    inf_induced_hosp_risk_reduce = params.inf_induced_hosp_risk_reduce
    inf_induced_proportional_risk_reduce = inf_induced_hosp_risk_reduce / (1 - inf_induced_hosp_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_hosp_risk_reduce_initial

    prob_hosp = (params.IP_to_ISH_prop / immunity_force) * vax_immunity_factor

    rate = params.IP_to_IS_rate * (1 - prob_hosp)

    return rate

compute_ISH_to_HD_rate(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_ISH_to_HD_rate(state: FluFullMetapopStateTensors,
                           params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    inf_induced_death_risk_reduce = params.inf_induced_death_risk_reduce

    inf_induced_proportional_risk_reduce = \
        inf_induced_death_risk_reduce / (1 - inf_induced_death_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_death_risk_reduce_initial

    prob_death = (params.ISH_to_HD_prop / immunity_force) * vax_immunity_factor

    rate = prob_death * params.ISH_to_H_rate

    return rate

compute_ISH_to_HR_rate(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_ISH_to_HR_rate(state: FluFullMetapopStateTensors,
                           params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    inf_induced_death_risk_reduce = params.inf_induced_death_risk_reduce

    inf_induced_proportional_risk_reduce = \
        inf_induced_death_risk_reduce / (1 - inf_induced_death_risk_reduce)

    immunity_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_death_risk_reduce_initial

    prob_death = (params.ISH_to_HD_prop / immunity_force) * vax_immunity_factor

    rate = (1 - prob_death) * params.ISH_to_H_rate

    return rate

compute_ISR_to_R(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_ISR_to_R(state: FluFullMetapopStateTensors,
                     params: FluFullMetapopParamsTensors,
                     dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    rate = params.ISR_to_R_rate

    ISR_to_R = state.ISR * torch_approx_binom_probability_from_rate(rate, dt)

    return ISR_to_R

compute_MV_change(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, precomputed: FluPrecomputedTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_MV_change(state: FluFullMetapopStateTensors,
                      params: FluFullMetapopParamsTensors,
                      precomputed: FluPrecomputedTensors,
                      dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    MV_change = state.daily_vaccines - \
                params.vax_induced_immune_wane * state.MV

    return MV_change * dt

compute_M_change(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, precomputed: FluPrecomputedTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_M_change(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors,
                     precomputed: FluPrecomputedTensors,
                     dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    # Note: already includes dt
    R_to_S = state.R * torch_approx_binom_probability_from_rate(params.R_to_S_rate, dt)

    M_change = (R_to_S / precomputed.total_pop_LAR_tensor) * \
               (1 - params.inf_induced_saturation * state.M) - \
               params.inf_induced_immune_wane * state.M * dt

    # Because R_to_S includes dt already, we do not return M_change * dt -- we only multiply
    #   the last term in the expression above by dt
    return M_change

compute_R_to_S(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, dt: float) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_R_to_S(state: FluFullMetapopStateTensors,
                   params: FluFullMetapopParamsTensors,
                   dt: float) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))
    """

    rate = params.R_to_S_rate

    R_to_S = state.R * torch_approx_binom_probability_from_rate(rate, dt)

    return R_to_S

compute_S_to_E(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, precomputed: FluPrecomputedTensors, dt: float, total_mixing_exposure: torch.Tensor = None) -> torch.Tensor

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

If total_mixing_exposure is provided, use it directly (for daily-update mode matching the numpy metapop model). Otherwise compute it from current state.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_S_to_E(state: FluFullMetapopStateTensors,
                   params: FluFullMetapopParamsTensors,
                   precomputed: FluPrecomputedTensors,
                   dt: float,
                   total_mixing_exposure: torch.Tensor = None) -> torch.Tensor:
    """
    Returns:
        (torch.Tensor of size (L, A, R))

    If total_mixing_exposure is provided, use it directly (for daily-update mode
    matching the numpy metapop model). Otherwise compute it from current state.
    """

    if total_mixing_exposure is None:
        total_mixing_exposure = compute_total_mixing_exposure(state, params, precomputed)

    if total_mixing_exposure.size() != torch.Size([precomputed.L,
                                                   precomputed.A,
                                                   precomputed.R]):
        raise Exception("force_of_infection must be L x A x R corresponding \n"
                        "to number of locations (subpopulations), age groups, \n"
                        "and risk groups.")

    beta_adjusted = compute_beta_adjusted(state, params)

    inf_induced_inf_risk_reduce = params.inf_induced_inf_risk_reduce
    inf_induced_proportional_risk_reduce = inf_induced_inf_risk_reduce / (1 - inf_induced_inf_risk_reduce)

    immune_force = 1 + inf_induced_proportional_risk_reduce * state.M

    vax_immunity_factor = 1 - state.MV * params.vax_induced_inf_risk_reduce_initial

    rate = beta_adjusted * total_mixing_exposure * vax_immunity_factor / immune_force

    S_to_E = state.S * torch_approx_binom_probability_from_rate(rate, dt)

    return S_to_E

compute_active_pop_LAR(state: FluTravelStateTensors, _params: FluTravelParamsTensors, precomputed: FluPrecomputedTensors) -> torch.Tensor

Compute the active population for location-age-risk (l, a, r) as a tensor. Used to compute the effective population in the travel model, which is the population size adjusted for incoming visitors, residents traveling, and assuming hospitalized individuals are not mobile enough to infect others.

Returns:

Type Description
Tensor

torch.Tensor of size (L, A, R): Active population: those who are not hospitalized (i.e. those who are not too sick to move and travel regularly)

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_active_pop_LAR(state: FluTravelStateTensors,
                           _params: FluTravelParamsTensors,
                           precomputed: FluPrecomputedTensors) -> torch.Tensor:
    """
    Compute the active population for location-age-risk
    (l, a, r) as a tensor. Used to compute the
    effective population in the travel model, which is
    the population size adjusted for incoming visitors,
    residents traveling, and assuming hospitalized
    individuals are not mobile enough to infect others.

    Returns:
        torch.Tensor of size (L, A, R):
            Active population: those who are not
            hospitalized (i.e. those who are not too sick
            to move and travel regularly)
    """

    # _params is not used now -- but this is included for
    #   function signature consistency with other
    #   similar computation functions

    return precomputed.total_pop_LAR_tensor - state.HR - state.HD

compute_beta_adjusted(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors) -> torch.Tensor

Computes beta-adjusted humidity.

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, R))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_beta_adjusted(state: FluFullMetapopStateTensors,
                          params: FluFullMetapopParamsTensors) -> torch.Tensor:
    """
    Computes beta-adjusted humidity.

    Returns:
        (torch.Tensor of size (L, A, R))
    """

    absolute_humidity = state.absolute_humidity
    beta_adjusted = params.beta_baseline * (1 + params.humidity_impact * np.exp(-180 * absolute_humidity))

    return beta_adjusted

compute_effective_pop_LA(state: FluTravelStateTensors, params: FluTravelParamsTensors, precomputed: FluPrecomputedTensors) -> torch.Tensor

Returns:

Type Description
Tensor

torch.Tensor of size (L, A): Effective population, summed over risk groups. See compute_active_pop_LAR docstring for more information.

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_effective_pop_LA(state: FluTravelStateTensors,
                             params: FluTravelParamsTensors,
                             precomputed: FluPrecomputedTensors) -> torch.Tensor:
    """
    Returns:
        torch.Tensor of size (L, A):
            Effective population, summed over risk groups.
            See `compute_active_pop_LAR` docstring for more
            information.
    """

    active_pop_LAR = compute_active_pop_LAR(state, params, precomputed)

    # Nonlocal travel proportions is L x L
    # Active population LAR is L x A x R
    outside_visitors_LAR = torch.einsum("kl,kar->lar",
                                        precomputed.nonlocal_travel_prop,
                                        active_pop_LAR)

    # This is correct -- Dave checked in meeting -- we don't need Einstein
    #   notation here!
    # In computation, broadcast sum_residents_nonlocal_travel_prop to be L x 1 x 1
    traveling_residents_LAR = precomputed.sum_residents_nonlocal_travel_prop[:, None, None] * \
                              active_pop_LAR

    mobility_modifier = state.mobility_modifier[:, :, 0]

    effective_pop_LA = precomputed.total_pop_LA + mobility_modifier * \
                       torch.sum(outside_visitors_LAR - traveling_residents_LAR, dim=2)

    return effective_pop_LA

compute_flu_contact_matrix(params: FluFullMetapopParamsTensors, schedules: FluFullMetapopScheduleTensors, day_counter: int) -> torch.Tensor

Computes flu model contact matrix in tensor format -- makes adjustments based on whether day is school day or work day.

Returns:

Type Description
Tensor

(torch.Tensor of size (L, A, A))

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def compute_flu_contact_matrix(params: FluFullMetapopParamsTensors,
                               schedules: FluFullMetapopScheduleTensors,
                               day_counter: int) -> torch.Tensor:
    """
    Computes flu model contact matrix in tensor format -- makes
    adjustments based on whether day is school day or work day.

    Returns:
        (torch.Tensor of size (L, A, A))
    """

    # Here, using schedules.is_school_day[day_counter][:,:,0] and similarly for
    #   is_work_day because each contact matrix (as a metapop tensor) is L x A x A --
    #   we don't use risk -- assume here that we do not have a different school/work-day
    #   schedule based on risk, so just grab the first risk group
    # But then we have to take (1 - schedules.is_school_day[day_counter][:, :, 0]), which is
    #   L x A, and then make it L x A x 1 (unsqueeze the last dimension) to make the
    #   broadcasting work (because this gets element-wise multiplied by params.school_contact_matrix)
    flu_contact_matrix = \
        params.total_contact_matrix - \
        params.school_contact_matrix * (1 - schedules.is_school_day[day_counter][:, :, 0]).unsqueeze(dim=2) - \
        params.work_contact_matrix * (1 - schedules.is_work_day[day_counter][:, :, 0]).unsqueeze(dim=2)

    return flu_contact_matrix

compute_local_to_local_exposure(flu_contact_matrix: torch.Tensor, mobility_modifier: torch.Tensor, sum_residents_nonlocal_travel_prop: torch.Tensor, wtd_infectious_ratio_LLA: torch.Tensor, location_ix: int) -> torch.Tensor

Raw means that this is unnormalized by relative_suscept. Excludes beta and population-level immunity adjustments -- those are factored in later.

Returns:

Type Description
Tensor

torch.Tensor of size (A): For a given location (specified by location_ix), compute local transmission caused by residents traveling within their home location, summed over risk groups.

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_local_to_local_exposure(flu_contact_matrix: torch.Tensor,
                                    mobility_modifier: torch.Tensor,
                                    sum_residents_nonlocal_travel_prop: torch.Tensor,
                                    wtd_infectious_ratio_LLA: torch.Tensor,
                                    location_ix: int) -> torch.Tensor:
    """
    Raw means that this is unnormalized by `relative_suscept`.
    Excludes beta and population-level immunity adjustments --
    those are factored in later.

    Returns:
        torch.Tensor of size (A):
            For a given location (specified by `location_ix`), compute
            local transmission caused by residents traveling within their
            home location, summed over risk groups.
    """

    # WARNING: we assume `mobility_modifier` is input as (A, 1) --
    # if this changes, we have to change the implementation.
    # The risk dimension does not have unique values, so we just
    # grab the first element of the risk dimension.
    proportion_staying_home = np.maximum(
        0, 
        (1 - mobility_modifier[location_ix, :, 0] * sum_residents_nonlocal_travel_prop[location_ix])
        )

    result = torch.mul(
        proportion_staying_home,
        torch.matmul(
            flu_contact_matrix[location_ix, :, :],
            torch.mul(
                proportion_staying_home, 
                wtd_infectious_ratio_LLA[location_ix, location_ix, :]
                )
            )
        )

    return result

compute_outside_visitors_exposure(flu_contact_matrix: torch.Tensor, mobility_modifier: torch.Tensor, sum_residents_nonlocal_travel_prop: torch.Tensor, travel_proportions: torch.Tensor, wtd_infectious_ratio_LLA: torch.Tensor, local_ix: int, visitors_ix: int) -> torch.Tensor

Computes raw (unnormalized by relative_suscept) transmission to local_ix due to outside visitors from visitors_ix. Excludes beta and population-level immunity adjustments -- those are factored in later.

Returns:

Type Description
Tensor

torch.Tensor of size (A)

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_outside_visitors_exposure(flu_contact_matrix: torch.Tensor,
                                      mobility_modifier: torch.Tensor,
                                      sum_residents_nonlocal_travel_prop: torch.Tensor,
                                      travel_proportions: torch.Tensor,
                                      wtd_infectious_ratio_LLA: torch.Tensor,
                                      local_ix: int,
                                      visitors_ix: int) -> torch.Tensor:
    """
    Computes raw (unnormalized by `relative_suscept`) transmission
    to `local_ix` due to outside visitors from `visitors_ix`.
    Excludes beta and population-level immunity adjustments --
    those are factored in later.

    Returns:
        torch.Tensor of size (A)
    """

    # In location `local_ix`, we are looking at the visitors from
    #   `visitors_ix` who come to `local_ix` (and infect folks in `local_ix`)

    # See WARNING in `compute_local_to_local_exposure()`
    proportion_staying_home = np.maximum(
        0,
        (1 - mobility_modifier[local_ix, :, 0] * sum_residents_nonlocal_travel_prop[local_ix])
        )

    result = torch.mul(
        proportion_staying_home * travel_proportions[visitors_ix, local_ix], \
        torch.matmul(
            mobility_modifier[visitors_ix, :, 0] * flu_contact_matrix[local_ix, :, :],
            wtd_infectious_ratio_LLA[visitors_ix, local_ix, :]
            )
        )

    return result

compute_pop_by_age(subpop_params: FluSubpopParams) -> np.ndarray

Returns:

Type Description
ndarray

np.ndarray: A x 1 array -- where A is the number of age groups -- where ith element corresponds to total population (across all compartments, including "D", and across all risk groups) in age group i

Source code in CLT_BaseModel/flu_core/flu_components.py
def compute_pop_by_age(subpop_params: FluSubpopParams) -> np.ndarray:
    """
    Returns:
        np.ndarray:
            A x 1 array -- where A is the number of age groups --
            where ith element corresponds to total population
            (across all compartments, including "D", and across all risk groups)
            in age group i
    """

    return np.sum(subpop_params.total_pop_age_risk, axis=1, keepdims=True)

compute_residents_traveling_exposure(flu_contact_matrix: torch.Tensor, mobility_modifier: torch.Tensor, sum_residents_nonlocal_travel_prop: torch.Tensor, travel_proportions: torch.Tensor, wtd_infectious_ratio_LLA: torch.Tensor, local_ix: int, dest_ix: int) -> torch.Tensor

Computes raw (unnormalized by relative_suscept) transmission to local_ix, due to residents of local_ix traveling to dest_ix and getting infected in dest_ix. Excludes beta and population-level immunity adjustments -- those are factored in later.

Returns:

Type Description
Tensor

torch.Tensor of size (A)

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_residents_traveling_exposure(flu_contact_matrix: torch.Tensor,
                                         mobility_modifier: torch.Tensor,
                                         sum_residents_nonlocal_travel_prop: torch.Tensor,
                                         travel_proportions: torch.Tensor,
                                         wtd_infectious_ratio_LLA: torch.Tensor,
                                         local_ix: int,
                                         dest_ix: int) -> torch.Tensor:
    """
    Computes raw (unnormalized by `relative_suscept`) transmission
    to `local_ix`, due to residents of `local_ix` traveling to `dest_ix`
    and getting infected in `dest_ix`. Excludes beta and population-level
    immunity adjustments -- those are factored in later.

    Returns:
        torch.Tensor of size (A)
    """

    # See WARNING in `compute_local_to_local_exposure()`
    proportion_staying_home_at_dest = np.maximum(
        0,
        (1 - mobility_modifier[dest_ix, :, 0] * sum_residents_nonlocal_travel_prop[dest_ix])
        )

    # Vectorized version:
    # For k != dest_ix: infectious_proportion[k, a'] = mobility_modifier[k, a'] * travel_proportions[k, dest]
    # For k == dest_ix: infectious_proportion[k, a'] = proportion_staying_home_at_dest[a']
    # Then the sum collapses to: einsum("ka,ka->a", Q[:, dest, :], infectious_proportion)
    infectious_proportion = mobility_modifier[:, :, 0] * travel_proportions[:, dest_ix].unsqueeze(1)
    infectious_proportion[dest_ix, :] = proportion_staying_home_at_dest

    result = torch.mul(
        mobility_modifier[local_ix, :, 0] * travel_proportions[local_ix, dest_ix],
        torch.matmul(
            flu_contact_matrix[dest_ix, :, :],
            torch.einsum("ka,ka->a", wtd_infectious_ratio_LLA[:, dest_ix, :], infectious_proportion)
            )
        )

    return result

compute_rsquared(reference_timeseries: list[np.ndarray], simulated_timeseries: list[np.ndarray]) -> float

Source code in CLT_BaseModel/flu_core/flu_accept_reject.py
def compute_rsquared(reference_timeseries: list[np.ndarray],
                     simulated_timeseries: list[np.ndarray]) -> float:
    if len(reference_timeseries) != len(simulated_timeseries):
        raise ValueError("Reference time series and simulated time series \n"
                         "must have same length.")

    reference_timeseries = np.asarray(reference_timeseries)
    simulated_timeseries = np.asarray(simulated_timeseries)

    ybar = reference_timeseries.mean(axis=0)

    ss_residual = np.sum(np.square(simulated_timeseries - reference_timeseries))
    ss_total = np.sum(np.square(reference_timeseries - ybar))

    return 1 - ss_residual / ss_total

compute_total_mixing_exposure(state: FluTravelStateTensors, params: FluTravelParamsTensors, precomputed: FluPrecomputedTensors) -> torch.Tensor

Computes "total mixing exposure" for location-age-risk (l, a, r) -- the rate of exposure to infectious individuals, accounting for both local transmission, incoming visitors, and residents traveling. Normalized by relative_suscept!

Combines subroutines compute_local_to_local_exposure(), compute_outside_visitors_exposure(), and compute_residents_traveling_exposure(). Note that these subroutines do not include relative susceptibility -- but this function includes relative susceptibility -- this is to avoid unnecessary repeated multiplication by relative susceptible in each subroutine.

Returns:

Type Description
Tensor

torch.Tensor of size (L, A, R)

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_total_mixing_exposure(state: FluTravelStateTensors,
                                  params: FluTravelParamsTensors,
                                  precomputed: FluPrecomputedTensors) -> torch.Tensor:
    """
    Computes "total mixing exposure" for location-age-risk
    (l, a, r) -- the rate of exposure to infectious individuals,
    accounting for both local transmission, incoming visitors, and
    residents traveling. **Normalized by `relative_suscept`!**

    Combines subroutines `compute_local_to_local_exposure()`,
    `compute_outside_visitors_exposure()`, and `compute_residents_traveling_exposure()`.
    Note that these subroutines do not include relative susceptibility --
    but this function includes relative susceptibility -- this is to avoid
    unnecessary repeated multiplication by relative susceptible in each subroutine.

    Returns:
        torch.Tensor of size (L, A, R)
    """

    L, A, R = precomputed.L, precomputed.A, precomputed.R

    mobility_modifier = state.mobility_modifier
    flu_contact_matrix = state.flu_contact_matrix
    travel_proportions = params.travel_proportions

    sum_residents_nonlocal_travel_prop = precomputed.sum_residents_nonlocal_travel_prop
    wtd_infectious_ratio_LLA = compute_wtd_infectious_ratio_LLA(state, params, precomputed)

    relative_suscept = params.relative_suscept[0, :, 0]

    total_mixing_exposure = torch.tensor(np.zeros((L, A, R)))

    # Couldn't figure out how to do this without two for-loops ;)
    # Welcoming any efficiency improvements!
    for l in np.arange(L):

        raw_total_mixing_exposure = torch.tensor(np.zeros(A))

        raw_total_mixing_exposure = raw_total_mixing_exposure + \
                                    compute_local_to_local_exposure(flu_contact_matrix,
                                                                    mobility_modifier,
                                                                    sum_residents_nonlocal_travel_prop,
                                                                    wtd_infectious_ratio_LLA,
                                                                    l)

        for k in np.arange(L):
            if k == l:
                continue # no visit terms from a location to itself

            raw_total_mixing_exposure = raw_total_mixing_exposure + \
                                        compute_outside_visitors_exposure(
                                            flu_contact_matrix,
                                            mobility_modifier,
                                            sum_residents_nonlocal_travel_prop,
                                            travel_proportions,
                                            wtd_infectious_ratio_LLA,
                                            l,
                                            k)

            raw_total_mixing_exposure = raw_total_mixing_exposure + \
                                        compute_residents_traveling_exposure(
                                            flu_contact_matrix,
                                            mobility_modifier,
                                            sum_residents_nonlocal_travel_prop,
                                            travel_proportions,
                                            wtd_infectious_ratio_LLA,
                                            l,
                                            k)

        normalized_total_mixing_exposure = relative_suscept * raw_total_mixing_exposure

        total_mixing_exposure[l, :, :] = normalized_total_mixing_exposure.view(A, 1).expand((A, R))

    return total_mixing_exposure

compute_vax_induced_risk_reduce_initial(params: FluSubpopParams, schedules: sc.objdict, start_real_date: datetime.date) -> tuple

Computes the "zero-waning" (peak, just-after-protection-delay) vaccine efficacy values vax_induced_inf_risk_reduce_initial, vax_induced_hosp_risk_reduce_initial, and vax_induced_death_risk_reduce_initial from the corresponding season-average efficacy values (vax_induced_inf_risk_reduce, vax_induced_hosp_risk_reduce, vax_induced_death_risk_reduce), vax_induced_immune_wane, and the daily_vaccines schedule.

For each age-risk group, given the season's vaccination timing (p_prot, the proportion of effective doses -- i.e. vaccination date plus protection delay -- given on each day of the season) and waning rate w_V, this solves for the peak efficacy VE_0 such that the dose-timing-weighted average realized efficacy over the season equals the input (season-average) efficacy value:

VE_0 = w_V * T * VE_season /
       ((1 - exp(-w_V)) * sum_{tau=t0}^{T-1} [
           (sum_{u=t0}^{tau} p_prot(u) * exp(-w_V * (tau - u))) /
           (sum_{u=t0}^{tau} p_prot(u))
       ])

where t0 is the first day of the season with a nonzero dose and T is the number of days from t0 through the last day with a nonzero dose.

The season window is the period between two consecutive occurrences of vax_immunity_reset_date_mm_dd (the occurrence on or before start_real_date, through one year later), or if that parameter is not set, the 12 months starting from the first date covered by the daily_vaccines schedule. This window is further intersected with the actual date range covered by the daily_vaccines schedule (which matters when the reset-date window extends past the end of the schedule, or -- in the no-reset-date case -- when the schedule itself spans less than 12 months).

Edge cases
  • If waning (w_V) is 0 for a given age-risk group, VE_0 equals the input (season-average) value for that group -- no adjustment is needed since there is no waning to correct for.
  • If there are no vaccine doses in the (intersected) season window for a given age-risk group, VE_0 equals the input value for that group -- vaccine-induced immunity is always 0 for that group, so the value is never actually applied.
  • If params.adjust_VE_for_seasonal_waning is False, this adjustment is skipped entirely and the input (season-average) values are returned unchanged (broadcast to shape (A, R)).

Parameters:

Name Type Description Default
params FluSubpopParams

holds vax_induced_inf_risk_reduce, vax_induced_hosp_risk_reduce, vax_induced_death_risk_reduce, vax_induced_immune_wane, vax_immunity_reset_date_mm_dd, num_age_groups, and num_risk_groups.

required
schedules objdict

holds the daily_vaccines Schedule instance (already shifted by vax_protection_delay_days, date-indexed, one A x R array per day).

required
start_real_date date

real-world date corresponding to the start of the simulation -- used to anchor the season window.

required

Returns:

Type Description
tuple

(vax_induced_inf_risk_reduce_initial, vax_induced_hosp_risk_reduce_initial, vax_induced_death_risk_reduce_initial), each an np.ndarray of shape (A, R).

Source code in CLT_BaseModel/flu_core/flu_components.py
def compute_vax_induced_risk_reduce_initial(params: FluSubpopParams,
                                         schedules: sc.objdict,
                                         start_real_date: datetime.date) -> tuple:
    """
    Computes the "zero-waning" (peak, just-after-protection-delay) vaccine
    efficacy values `vax_induced_inf_risk_reduce_initial`,
    `vax_induced_hosp_risk_reduce_initial`, and `vax_induced_death_risk_reduce_initial`
    from the corresponding season-average efficacy values
    (`vax_induced_inf_risk_reduce`, `vax_induced_hosp_risk_reduce`,
    `vax_induced_death_risk_reduce`), `vax_induced_immune_wane`, and the
    `daily_vaccines` schedule.

    For each age-risk group, given the season's vaccination timing
    (`p_prot`, the proportion of effective doses -- i.e. vaccination date
    plus protection delay -- given on each day of the season) and waning
    rate w_V, this solves for the peak efficacy VE_0 such that the
    dose-timing-weighted average realized efficacy over the season equals
    the input (season-average) efficacy value:

        VE_0 = w_V * T * VE_season /
               ((1 - exp(-w_V)) * sum_{tau=t0}^{T-1} [
                   (sum_{u=t0}^{tau} p_prot(u) * exp(-w_V * (tau - u))) /
                   (sum_{u=t0}^{tau} p_prot(u))
               ])

    where t0 is the first day of the season with a nonzero dose and T is
    the number of days from t0 through the last day with a nonzero dose.

    The season window is the period between two consecutive occurrences of
    `vax_immunity_reset_date_mm_dd` (the occurrence on or before
    `start_real_date`, through one year later), or if that parameter is
    not set, the 12 months starting from the first date covered by the
    `daily_vaccines` schedule. This window is further intersected with
    the actual date range covered by the `daily_vaccines` schedule (which
    matters when the reset-date window extends past the end of the
    schedule, or -- in the no-reset-date case -- when the schedule itself
    spans less than 12 months).

    Edge cases:
        - If waning (w_V) is 0 for a given age-risk group, VE_0 equals the
          input (season-average) value for that group -- no adjustment is
          needed since there is no waning to correct for.
        - If there are no vaccine doses in the (intersected) season window
          for a given age-risk group, VE_0 equals the input value for that
          group -- vaccine-induced immunity is always 0 for that group, so
          the value is never actually applied.
        - If `params.adjust_VE_for_seasonal_waning` is False, this
          adjustment is skipped entirely and the input (season-average)
          values are returned unchanged (broadcast to shape (A, R)).

    Args:
        params (FluSubpopParams):
            holds `vax_induced_inf_risk_reduce`, `vax_induced_hosp_risk_reduce`,
            `vax_induced_death_risk_reduce`, `vax_induced_immune_wane`,
            `vax_immunity_reset_date_mm_dd`, `num_age_groups`, and
            `num_risk_groups`.
        schedules (sc.objdict):
            holds the `daily_vaccines` `Schedule` instance (already shifted
            by `vax_protection_delay_days`, date-indexed, one A x R array
            per day).
        start_real_date (datetime.date):
            real-world date corresponding to the start of the simulation --
            used to anchor the season window.

    Returns:
        (vax_induced_inf_risk_reduce_initial, vax_induced_hosp_risk_reduce_initial,
         vax_induced_death_risk_reduce_initial), each an np.ndarray of shape (A, R).
    """

    target_shape = (params.num_age_groups, params.num_risk_groups)

    if params.adjust_VE_for_seasonal_waning is False:
        return (
            np.broadcast_to(np.asarray(params.vax_induced_inf_risk_reduce, dtype=float), target_shape).copy(),
            np.broadcast_to(np.asarray(params.vax_induced_hosp_risk_reduce, dtype=float), target_shape).copy(),
            np.broadcast_to(np.asarray(params.vax_induced_death_risk_reduce, dtype=float), target_shape).copy(),
        )

    w_arr = np.broadcast_to(
        np.asarray(params.vax_induced_immune_wane, dtype=float), target_shape)

    vaccines_df = schedules["daily_vaccines"].timeseries_df

    schedule_min_date = vaccines_df.index.min()
    schedule_max_date = vaccines_df.index.max()

    if params.vax_immunity_reset_date_mm_dd is not None:
        month, day = (int(x) for x in params.vax_immunity_reset_date_mm_dd.split('_'))
        window_start = datetime.date(start_real_date.year, month, day)
        if window_start >= start_real_date:
            window_start = datetime.date(start_real_date.year - 1, month, day)
        window_end = datetime.date(window_start.year + 1, month, day)
    else:
        window_start = schedule_min_date
        window_end = schedule_min_date + datetime.timedelta(days=365)

    window_start = max(window_start, schedule_min_date)
    window_end = min(window_end, schedule_max_date + datetime.timedelta(days=1))

    if window_start < window_end:
        mask = (vaccines_df.index >= window_start) & (vaccines_df.index < window_end)
        window_doses_df = vaccines_df.loc[mask]
    else:
        window_doses_df = vaccines_df.iloc[0:0]

    if not window_doses_df.empty:
        doses_stack = np.stack(window_doses_df["daily_vaccines"].values, axis=0)
        if doses_stack.shape[1:] != target_shape:
            # Time series has a different age-risk resolution than the
            # risk-reduce parameters -- aggregate (sum) across all
            # dimensions and broadcast the resulting total evenly across
            # every age-risk group.
            doses_stack = doses_stack.reshape(doses_stack.shape[0], -1).sum(axis=1, keepdims=True)
            doses_stack = np.broadcast_to(doses_stack, (doses_stack.shape[0],) + target_shape)
    else:
        doses_stack = np.zeros((0,) + target_shape)

    def compute_initial_for_field(field_name: str) -> np.ndarray:

        ve_season_arr = np.broadcast_to(
            np.asarray(getattr(params, field_name), dtype=float), target_shape).copy()
        ve_initial_arr = ve_season_arr.copy()

        for a in range(target_shape[0]):
            for r in range(target_shape[1]):

                w_ar = w_arr[a, r]
                ve_season = ve_season_arr[a, r]

                if w_ar == 0 or doses_stack.shape[0] == 0:
                    continue

                cell_doses = doses_stack[:, a, r]
                nonzero_idxs = np.flatnonzero(cell_doses > 0)

                if nonzero_idxs.size == 0:
                    continue

                t0_idx, last_idx = nonzero_idxs[0], nonzero_idxs[-1]
                sub = cell_doses[t0_idx:last_idx + 1]
                T = sub.size
                p_prot = sub / sub.sum()

                cumsum_p_prot = np.cumsum(p_prot)

                # numer[n] = sum_{u=0}^{n} p_prot(u) * exp(-w * (n - u))
                #          = p_prot(n) + exp(-w) * numer[n - 1]
                decay = np.exp(-w_ar)
                numer = np.empty(T)
                acc = 0.0
                for n in range(T):
                    acc = p_prot[n] + decay * acc
                    numer[n] = acc

                S = np.sum(numer / cumsum_p_prot)

                ve_initial_arr[a, r] = w_ar * T * ve_season / ((1 - decay) * S)

        return ve_initial_arr

    return (compute_initial_for_field("vax_induced_inf_risk_reduce"),
           compute_initial_for_field("vax_induced_hosp_risk_reduce"),
           compute_initial_for_field("vax_induced_death_risk_reduce"))

compute_wtd_infectious_LA(state: FluTravelStateTensors, params: FluTravelParamsTensors) -> torch.Tensor

Returns:

Type Description
Tensor

torch.Tensor of size (L, A): Weighted infectious, summed over risk groups: includes presymptomatic, asymptomatic, and symptomatic, weighted by relative infectiousness

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_wtd_infectious_LA(state: FluTravelStateTensors,
                              params: FluTravelParamsTensors) -> torch.Tensor:
    """
    Returns:
        torch.Tensor of size (L, A):
            Weighted infectious, summed over risk groups:
            includes presymptomatic, asymptomatic, and symptomatic,
            weighted by relative infectiousness
    """

    # Einstein notation here means sum over risk groups
    ISR = torch.einsum("lar->la", state.ISR)
    ISH = torch.einsum("lar->la", state.ISH)
    wtd_IP = \
        params.IP_relative_inf * torch.einsum("lar->la", state.IP)
    wtd_IA = \
        params.IA_relative_inf * torch.einsum("lar->la", state.IA)

    return ISR + ISH + wtd_IP + wtd_IA

compute_wtd_infectious_ratio_LLA(state: FluTravelStateTensors, params: FluTravelParamsTensors, precomputed: FluPrecomputedTensors) -> torch.Tensor

Returns:

Type Description
Tensor

torch.Tensor of size (L, L, A): Element i,j,a corresponds to ratio of weighted infectious people in location i, age group a (summed over risk groups) to the effective population in location j (summed over risk groups)

Source code in CLT_BaseModel/flu_core/flu_travel_functions.py
def compute_wtd_infectious_ratio_LLA(state: FluTravelStateTensors,
                                     params: FluTravelParamsTensors,
                                     precomputed: FluPrecomputedTensors) -> torch.Tensor:
    """
    Returns:
        torch.Tensor of size (L, L, A):
            Element i,j,a corresponds to ratio of weighted infectious people
            in location i, age group a (summed over risk groups) to the effective
            population in location j (summed over risk groups)
    """

    wtd_infectious_LA = compute_wtd_infectious_LA(state, params)

    effective_pop_LA = compute_effective_pop_LA(state, params, precomputed)

    prop_wtd_infectious = torch.einsum("ka,la->kla",
                                       wtd_infectious_LA,
                                       1 / effective_pop_LA)

    return prop_wtd_infectious

compute_wtd_presymp_asymp_by_age(subpop_state: FluSubpopState, subpop_params: FluSubpopParams) -> np.ndarray

Returns weighted sum of IP and IA compartment for subpopulation with given state and parameters. IP and IA are weighted by their relative infectiousness respectively, and then summed over risk groups.

Returns:

Type Description
ndarray

np.ndarray of shape (A, R)

Source code in CLT_BaseModel/flu_core/flu_components.py
def compute_wtd_presymp_asymp_by_age(subpop_state: FluSubpopState,
                                     subpop_params: FluSubpopParams) -> np.ndarray:
    """
    Returns weighted sum of IP and IA compartment for
        subpopulation with given state and parameters.
        IP and IA are weighted by their relative infectiousness
        respectively, and then summed over risk groups.

    Returns:
        np.ndarray of shape (A, R)
    """

    # sum over risk groups
    wtd_IP = \
        subpop_params.IP_relative_inf * np.sum(subpop_state.IP, axis=1, keepdims=True)
    wtd_IA = \
        subpop_params.IA_relative_inf * np.sum(subpop_state.IA, axis=1, keepdims=True)

    return wtd_IP + wtd_IA

create_dict_of_tensors(d: dict, requires_grad: bool = True) -> dict

Converts dictionary entries to tensor (of type torch.float32) and if requires_grad is True, turns on gradient tracking for each entry -- returns new dictionary.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def create_dict_of_tensors(d: dict,
                           requires_grad: bool = True) -> dict:
    """
    Converts dictionary entries to `tensor` (of type `torch.float32`)
    and if `requires_grad` is `True`, turns on gradient tracking for
    each entry -- returns new dictionary.
    """

    def to_tensor(k, v):
        if v is None:
            return None
        else:
            return torch.tensor(v, dtype=torch.float32, requires_grad=requires_grad)

    return {k: to_tensor(k, v) for k, v in d.items()}

create_timeseries_df_from_day_of_week_schedule(day_of_week_schedule: pd.DataFrame, start_date: datetime.date) -> pd.DataFrame

Creates a dataframe containing a timeseries of values for each date starting from start_date for 10 years.

Parameters

day_of_week_schedule : pd.DataFrame Column day_of_week with values monday, tuesday, ... Second column has values for that day of week. start_date : datetime.date First day in timeseries.

Returns

pd.DataFrame Column date with all dates from start_date for 10 years. Second column has values for that date.

Source code in CLT_BaseModel/flu_core/flu_components.py
def create_timeseries_df_from_day_of_week_schedule(
        day_of_week_schedule: pd.DataFrame,
        start_date: datetime.date) -> pd.DataFrame:
    """
    Creates a dataframe containing a timeseries of values
    for each date starting from start_date for 10 years.

    Parameters
    ----------
    day_of_week_schedule : pd.DataFrame
        Column day_of_week with values monday, tuesday, ...
        Second column has values for that day of week.
    start_date : datetime.date
        First day in timeseries.

    Returns
    -------
    pd.DataFrame
        Column date with all dates from start_date for 10 years.
        Second column has values for that date.
    """

    df_day_of_week = day_of_week_schedule.copy()

    # Create full timeseries dataframe by repeating day of week schedule
    duration_days = 10 * 365 # extend to 10 years to be safe
    new_dates = pd.date_range(start=start_date, periods=duration_days, freq='D')
    df = pd.DataFrame({'date': new_dates})

    df['day_of_week'] = df['date'].dt.day_name().str.lower()
    df = pd.merge(
        df, df_day_of_week, 
        on='day_of_week', how='left'
        ).drop(columns=['day_of_week'])

    df = df.set_index('date')

    return df

cumulative_deaths(metapop_model: MetapopModel, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None) -> float

Season-total deaths = sum of HD→D transition flows (scalar).

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def cumulative_deaths(
    metapop_model: MetapopModel,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
) -> float:
    """
    Season-total deaths = sum of HD→D transition flows (scalar).
    """
    arr = _tvar_daily(metapop_model, ["HD_to_D"], subpop_name)
    return float(_apply_ar_filter(arr, age_group, risk_group).sum())

cumulative_hospitalizations(metapop_model: MetapopModel, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None) -> float

Season-total hospital admissions (scalar).

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def cumulative_hospitalizations(
    metapop_model: MetapopModel,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
) -> float:
    """
    Season-total hospital admissions (scalar).
    """
    return float(
        daily_hospital_admissions(metapop_model, subpop_name, age_group, risk_group).sum()
    )

daily_deaths(metapop_model: MetapopModel, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None) -> np.ndarray

Daily deaths = HD→D transition flows, aggregated to daily totals.

Returns

np.ndarray, shape (days,)

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def daily_deaths(
    metapop_model: MetapopModel,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
) -> np.ndarray:
    """
    Daily deaths = HD→D transition flows, aggregated to daily totals.

    Returns
    -------
    np.ndarray, shape (days,)
    """
    arr = _tvar_daily(metapop_model, ["HD_to_D"], subpop_name)
    return _apply_ar_filter(arr, age_group, risk_group)

daily_hospital_admissions(metapop_model: MetapopModel, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None) -> np.ndarray

Daily new hospital admissions.

Computed as the sum of ISH→HR and ISH→HD transition flows, summed across subpopulations and timesteps-per-day.

Parameters

metapop_model : MetapopModel subpop_name : str or None Restrict to one subpopulation; None sums all. age_group : int or None Index of age group; None sums all age groups. risk_group : int or None Index of risk group; None sums all risk groups.

Returns

np.ndarray, shape (days,)

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def daily_hospital_admissions(
    metapop_model: MetapopModel,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
) -> np.ndarray:
    """
    Daily new hospital admissions.

    Computed as the sum of ISH→HR and ISH→HD transition flows,
    summed across subpopulations and timesteps-per-day.

    Parameters
    ----------
    metapop_model : MetapopModel
    subpop_name : str or None
        Restrict to one subpopulation; None sums all.
    age_group : int or None
        Index of age group; None sums all age groups.
    risk_group : int or None
        Index of risk group; None sums all risk groups.

    Returns
    -------
    np.ndarray, shape (days,)
    """
    arr = _tvar_daily(metapop_model, ["ISH_to_HR", "ISH_to_HD"], subpop_name)
    return _apply_ar_filter(arr, age_group, risk_group)

daily_new_infections(metapop_model: MetapopModel, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None) -> np.ndarray

Daily new infections = S→E transition flows, aggregated to daily totals.

Returns

np.ndarray, shape (days,)

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def daily_new_infections(
    metapop_model: MetapopModel,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
) -> np.ndarray:
    """
    Daily new infections = S→E transition flows, aggregated to daily totals.

    Returns
    -------
    np.ndarray, shape (days,)
    """
    arr = _tvar_daily(metapop_model, ["S_to_E"], subpop_name)
    return _apply_ar_filter(arr, age_group, risk_group)

plot_attack_rate_by_age(metapop_model: MetapopModel, ax: matplotlib.axes.Axes = None, savefig_filename: str = None, subpop_name: Optional[str] = None, title: str = None) -> matplotlib.axes.Axes

Bar chart of attack rate per age group.

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def plot_attack_rate_by_age(
    metapop_model: MetapopModel,
    ax: matplotlib.axes.Axes = None,
    savefig_filename: str = None,
    subpop_name: Optional[str] = None,
    title: str = None,
) -> matplotlib.axes.Axes:
    """
    Bar chart of attack rate per age group.
    """
    ax_provided = ax is not None
    if ax is None:
        fig, ax = plt.subplots(figsize=(8, 4))

    num_age_groups = _first_subpop(metapop_model).params.num_age_groups
    rates = [
        attack_rate(metapop_model, subpop_name, age_group=a)
        for a in range(num_age_groups)
    ]

    ax.bar(range(num_age_groups), rates, alpha=0.8)
    ax.set_xlabel("Age group")
    ax.set_ylabel("Attack rate")
    ax.set_xticks(range(num_age_groups))
    ax.set_xticklabels([f"Age {a}" for a in range(num_age_groups)])
    ax.set_title(title or "Attack rate by age group")

    if savefig_filename:
        plt.savefig(savefig_filename, dpi=300, bbox_inches="tight")
    if not ax_provided:
        plt.tight_layout()
        plt.show()

    return ax

plot_compartment_history(metapop_model: MetapopModel, compartment_names=('S', 'E', 'IP', 'ISR', 'ISH', 'IA', 'HR', 'HD', 'R', 'D'), ax: matplotlib.axes.Axes = None, savefig_filename: str = None, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None, title: str = None, linestyle: str = '-', label_suffix: str = '') -> matplotlib.axes.Axes

Time series of selected compartments, aggregated or stratified by age/risk.

Replaces the unsegmented plot_metapop_basic_compartment_history.

Parameters

metapop_model : MetapopModel compartment_names : sequence of str Compartments to include. ax : matplotlib.axes.Axes or None Axis to draw on; created if not provided. savefig_filename : str or None subpop_name : str or None Restrict to one subpopulation; None combines all. age_group : int or None When an integer, restrict to that specific age group. risk_group : int or None Filter to a single risk group (summed otherwise). title : str or None linestyle : str Matplotlib linestyle (e.g. "-", "--", ":"). label_suffix : str Appended to each legend entry (e.g. "[beta=0.05]").

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def plot_compartment_history(
    metapop_model: MetapopModel,
    compartment_names=("S", "E", "IP", "ISR", "ISH", "IA", "HR", "HD", "R", "D"),
    ax: matplotlib.axes.Axes = None,
    savefig_filename: str = None,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
    title: str = None,
    linestyle: str = "-",
    label_suffix: str = "",
) -> matplotlib.axes.Axes:
    """
    Time series of selected compartments, aggregated or stratified by age/risk.

    Replaces the unsegmented ``plot_metapop_basic_compartment_history``.

    Parameters
    ----------
    metapop_model : MetapopModel
    compartment_names : sequence of str
        Compartments to include.
    ax : matplotlib.axes.Axes or None
        Axis to draw on; created if not provided.
    savefig_filename : str or None
    subpop_name : str or None
        Restrict to one subpopulation; None combines all.
    age_group : int or None
        When an integer, restrict to that specific age group.
    risk_group : int or None
        Filter to a single risk group (summed otherwise).
    title : str or None
    linestyle : str
        Matplotlib linestyle (e.g. ``"-"``, ``"--"``, ``":"``).
    label_suffix : str
        Appended to each legend entry (e.g. `` "[beta=0.05]"``).
    """
    ax_provided = ax is not None
    if ax is None:
        fig, ax = plt.subplots(figsize=(10, 5))

    if subpop_name is not None:
        subpops = [metapop_model.subpop_models[subpop_name]]
    else:
        subpops = list(metapop_model.subpop_models.values())

    for comp_name in compartment_names:
        arrays = [
            np.asarray(sp.compartments[comp_name].history_vals_list)
            for sp in subpops
        ]
        total = np.sum(np.stack(arrays, axis=0), axis=0)  # (T, A, R)

        if age_group is not None:
            if risk_group is not None:
                vals = total[:, age_group, risk_group]
            else:
                vals = total[:, age_group, :].sum(axis=1)
        else:
            if risk_group is not None:
                vals = total[:, :, risk_group].sum(axis=1)
            else:
                vals = total.sum(axis=(1, 2))
        ax.plot(vals, label=f"{comp_name}{label_suffix}", alpha=0.7, linestyle=linestyle)

    ax.set_xlabel("Simulation day")
    ax.set_ylabel("Number of individuals")
    ax.set_title(title or "Compartment histories")
    ax.legend(fontsize=7)

    if savefig_filename:
        plt.savefig(savefig_filename, dpi=300, bbox_inches="tight")
    if not ax_provided:
        plt.tight_layout()
        plt.show()

    return ax

plot_daily_hospital_admissions(models, ax: matplotlib.axes.Axes = None, savefig_filename: str = None, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None, title: str = None) -> matplotlib.axes.Axes

Daily ISH→HR + ISH→HD, optionally overlaid across multiple models.

Parameters

models : MetapopModel or dict[str, MetapopModel] A single model, or {scenario_name: model} dict for overlay. For multi-replicate data, pass {scenario_name: list_of_MetapopModels} β€” in that case the median line and shaded 95 % interval are shown.

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def plot_daily_hospital_admissions(
    models,
    ax: matplotlib.axes.Axes = None,
    savefig_filename: str = None,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
    title: str = None,
) -> matplotlib.axes.Axes:
    """
    Daily ISH→HR + ISH→HD, optionally overlaid across multiple models.

    Parameters
    ----------
    models : MetapopModel  or  dict[str, MetapopModel]
        A single model, or ``{scenario_name: model}`` dict for overlay.
        For multi-replicate data, pass
        ``{scenario_name: list_of_MetapopModels}`` β€” in that case the median
        line and shaded 95 % interval are shown.
    """
    ax_provided = ax is not None
    if ax is None:
        fig, ax = plt.subplots(figsize=(10, 4))

    if isinstance(models, dict):
        prop_cycle = plt.rcParams["axes.prop_cycle"].by_key()["color"]
        for i, (scenario_name, model_or_list) in enumerate(models.items()):
            color = prop_cycle[i % len(prop_cycle)]
            if isinstance(model_or_list, list):
                # Multi-replicate: median + 95 % CI ribbon
                all_vals = np.stack(
                    [daily_hospital_admissions(m, subpop_name, age_group, risk_group)
                     for m in model_or_list],
                    axis=0,
                )  # (reps, days)
                median = np.median(all_vals, axis=0)
                lo = np.percentile(all_vals, 2.5, axis=0)
                hi = np.percentile(all_vals, 97.5, axis=0)
                days = np.arange(median.shape[0])
                ax.plot(days, median, label=scenario_name, color=color, alpha=0.9)
                ax.fill_between(days, lo, hi, color=color, alpha=0.2)
            else:
                vals = daily_hospital_admissions(
                    model_or_list, subpop_name, age_group, risk_group
                )
                ax.plot(vals, label=scenario_name, color=color, alpha=0.8)
    else:
        vals = daily_hospital_admissions(models, subpop_name, age_group, risk_group)
        ax.plot(vals, alpha=0.8)

    ax.set_xlabel("Day")
    ax.set_ylabel("Daily hospital admissions")
    ax.set_title(title or "Daily hospital admissions")
    ax.legend()

    if savefig_filename:
        plt.savefig(savefig_filename, dpi=300, bbox_inches="tight")
    if not ax_provided:
        plt.tight_layout()
        plt.show()

    return ax

plot_daily_new_infections(metapop_model: MetapopModel, ax: matplotlib.axes.Axes = None, savefig_filename: str = None, subpop_name: Optional[str] = None, age_group: Optional[int] = None, risk_group: Optional[int] = None, label: str = None, color=None, title: str = None) -> matplotlib.axes.Axes

Daily S→E flow aggregated to one curve per subpopulation (or combined).

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def plot_daily_new_infections(
    metapop_model: MetapopModel,
    ax: matplotlib.axes.Axes = None,
    savefig_filename: str = None,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    risk_group: Optional[int] = None,
    label: str = None,
    color=None,
    title: str = None,
) -> matplotlib.axes.Axes:
    """
    Daily S→E flow aggregated to one curve per subpopulation (or combined).
    """
    ax_provided = ax is not None
    if ax is None:
        fig, ax = plt.subplots(figsize=(10, 4))

    vals = daily_new_infections(metapop_model, subpop_name, age_group, risk_group)

    kw = {"alpha": 0.8}
    if color is not None:
        kw["color"] = color
    ax.plot(vals, label=label or "New infections", **kw)

    ax.set_xlabel("Day")
    ax.set_ylabel("Daily new infections")
    ax.set_title(title or "Daily new infections")
    ax.legend()

    if savefig_filename:
        plt.savefig(savefig_filename, dpi=300, bbox_inches="tight")
    if not ax_provided:
        plt.tight_layout()
        plt.show()

    return ax

plot_epi_metrics(metapop_model: MetapopModel, metric_names=('M', 'MV'), ax: matplotlib.axes.Axes = None, savefig_filename: str = None, subpop_name: Optional[str] = None, age_group: Optional[int] = None, title: str = None, linestyle: str = '-', label_suffix: str = '') -> matplotlib.axes.Axes

M and MV over time, restricted to age_group when provided.

Parameters

metapop_model : MetapopModel metric_names : sequence of str Epi metrics to plot; subset of ("M", "MV"). age_group : int or None When not None, restrict to that specific age group. linestyle : str Matplotlib linestyle (e.g. "-", "--", ":"). label_suffix : str Appended to each legend entry.

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def plot_epi_metrics(
    metapop_model: MetapopModel,
    metric_names=("M", "MV"),
    ax: matplotlib.axes.Axes = None,
    savefig_filename: str = None,
    subpop_name: Optional[str] = None,
    age_group: Optional[int] = None,
    title: str = None,
    linestyle: str = "-",
    label_suffix: str = "",
) -> matplotlib.axes.Axes:
    """
    M and MV over time, restricted to ``age_group`` when provided.

    Parameters
    ----------
    metapop_model : MetapopModel
    metric_names : sequence of str
        Epi metrics to plot; subset of ``("M", "MV")``.
    age_group : int or None
        When not None, restrict to that specific age group.
    linestyle : str
        Matplotlib linestyle (e.g. ``"-"``, ``"--"``, ``":"``).
    label_suffix : str
        Appended to each legend entry.
    """
    ax_provided = ax is not None
    if ax is None:
        fig, ax = plt.subplots(figsize=(10, 4))

    if subpop_name is not None:
        subpops = [metapop_model.subpop_models[subpop_name]]
    else:
        subpops = list(metapop_model.subpop_models.values())

    for metric_name in metric_names:
        arrays = [
            np.asarray(sp.epi_metrics[metric_name].history_vals_list)
            for sp in subpops
        ]
        # Average across subpops (epi metrics are fractions, not counts)
        total = np.mean(np.stack(arrays, axis=0), axis=0)  # (T, A, R)

        if age_group is not None:
            vals = total[:, age_group, :].mean(axis=1)
        else:
            vals = total.mean(axis=(1, 2))
        ax.plot(vals, label=f"{metric_name}{label_suffix}", alpha=0.7, linestyle=linestyle)

    ax.set_xlabel("Simulation day")
    ax.set_ylabel("Immunity level")
    ax.set_title(title or "Epi metrics (M, MV)")
    ax.legend()

    if savefig_filename:
        plt.savefig(savefig_filename, dpi=300, bbox_inches="tight")
    if not ax_provided:
        plt.tight_layout()
        plt.show()

    return ax

plot_scenario_comparison(models_dict: dict, metric_fn: Callable, ax: matplotlib.axes.Axes = None, savefig_filename: str = None, metric_name: str = None, title: str = None, **metric_kwargs) -> matplotlib.axes.Axes

Bar or box plot comparing a scalar metric across scenarios.

Parameters

models_dict : dict {scenario_name: model} for single-run bar chart, or {scenario_name: [model_rep1, model_rep2, ...]} for multi-replicate box plots. metric_fn : callable Scalar metric function (e.g. cumulative_hospitalizations). metric_name : str or None Label for the y-axis; defaults to metric_fn.__name__. **metric_kwargs Forwarded to metric_fn.

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def plot_scenario_comparison(
    models_dict: dict,
    metric_fn: Callable,
    ax: matplotlib.axes.Axes = None,
    savefig_filename: str = None,
    metric_name: str = None,
    title: str = None,
    **metric_kwargs,
) -> matplotlib.axes.Axes:
    """
    Bar or box plot comparing a scalar metric across scenarios.

    Parameters
    ----------
    models_dict : dict
        ``{scenario_name: model}`` for single-run bar chart, or
        ``{scenario_name: [model_rep1, model_rep2, ...]}`` for multi-replicate
        box plots.
    metric_fn : callable
        Scalar metric function (e.g. ``cumulative_hospitalizations``).
    metric_name : str or None
        Label for the y-axis; defaults to ``metric_fn.__name__``.
    **metric_kwargs
        Forwarded to ``metric_fn``.
    """
    ax_provided = ax is not None
    if ax is None:
        fig, ax = plt.subplots(figsize=(8, 5))

    names = list(models_dict.keys())
    values = []
    multi_rep = False

    for model_or_list in models_dict.values():
        if isinstance(model_or_list, list):
            multi_rep = True
            values.append([metric_fn(m, **metric_kwargs) for m in model_or_list])
        else:
            values.append(metric_fn(model_or_list, **metric_kwargs))

    if multi_rep:
        ax.boxplot(values, labels=names)
    else:
        ax.bar(names, values, alpha=0.8)

    ax.set_ylabel(metric_name or metric_fn.__name__)
    ax.set_title(title or "Scenario comparison")

    if savefig_filename:
        plt.savefig(savefig_filename, dpi=300, bbox_inches="tight")
    if not ax_provided:
        plt.tight_layout()
        plt.show()

    return ax

summarize_outcomes(values, credible_interval: float = 0.95) -> dict

Summarize a (reps,) array of scalar outcomes across replicates.

Parameters

values : array-like, shape (reps,) credible_interval : float Width of the central credible interval (default 0.95 β†’ 2.5th–97.5th percentiles).

Returns

dict with keys: mean, median, lower_ci, upper_ci

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def summarize_outcomes(
    values,
    credible_interval: float = 0.95,
) -> dict:
    """
    Summarize a (reps,) array of scalar outcomes across replicates.

    Parameters
    ----------
    values : array-like, shape (reps,)
    credible_interval : float
        Width of the central credible interval (default 0.95 β†’ 2.5th–97.5th
        percentiles).

    Returns
    -------
    dict with keys: ``mean``, ``median``, ``lower_ci``, ``upper_ci``
    """
    values = np.asarray(values, dtype=float)
    half = (1.0 - credible_interval) / 2.0
    return {
        "mean":     float(np.mean(values)),
        "median":   float(np.median(values)),
        "lower_ci": float(np.percentile(values, 100 * half)),
        "upper_ci": float(np.percentile(values, 100 * (1.0 - half))),
    }

torch_approx_binom_probability_from_rate(rate, dt)

Torch-compatible implementation of converting a rate into a probability. See analogous numpy implementation base_components/approx_binom_probability_from_rate() docstring for details.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def torch_approx_binom_probability_from_rate(rate, dt):
    """
    Torch-compatible implementation of converting a
    rate into a probability. See analogous numpy implementation
    `base_components/approx_binom_probability_from_rate()` docstring
    for details.
    """

    return 1 - torch.exp(-rate * dt)

torch_simulate_full_history(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, precomputed: FluPrecomputedTensors, schedules: FluFullMetapopScheduleTensors, num_days: int, timesteps_per_day: int) -> Tuple[dict, dict]

Simulates the flu model with a differentiable torch implementation that carries out binom_deterministic_no_round transition types -- returns hospital admits for calibration use.

See subroutine advance_timestep for additional details.

Returns:

Type Description
Tuple[dict, dict]

Returns compartment states and transition variables for day, location, age, risk, in tensor format.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def torch_simulate_full_history(state: FluFullMetapopStateTensors,
                                params: FluFullMetapopParamsTensors,
                                precomputed: FluPrecomputedTensors,
                                schedules: FluFullMetapopScheduleTensors,
                                num_days: int,
                                timesteps_per_day: int) -> Tuple[dict, dict]:
    """
    Simulates the flu model with a differentiable torch implementation
    that carries out `binom_deterministic_no_round` transition types --
    returns hospital admits for calibration use.

    See subroutine `advance_timestep` for additional details.

    Returns:
        (Tuple[dict, dict]):
            Returns compartment states and transition variables
            for day, location, age, risk, in tensor format.
    """

    dt = 1 / float(timesteps_per_day)

    state_history_dict = defaultdict(list)
    tvar_history_dict = defaultdict(list)

    # This could probably be written better so we don't have
    #   unused variables "_" that grab `advance_timestep` output?

    for day in range(num_days):
        state = update_state_with_schedules(state, params, schedules, day)
        # Compute mixing exposure once per day (matching numpy metapop model)
        daily_mixing_exposure = compute_total_mixing_exposure(state, params, precomputed)

        for timestep in range(timesteps_per_day):
            # TODO double check whether this split makes sense
            #   to get the total transition variables we should need to save values
            #   at each timestep when there are several steps per day
            #   (these variables may not be used anywhere right now)
            if timestep == timesteps_per_day-1:
                state, _, tvar_history = \
                    advance_timestep(state, params, precomputed, dt, save_tvar_history=True,
                                     total_mixing_exposure=daily_mixing_exposure)
                for key in tvar_history:
                    tvar_history_dict[key].append(tvar_history[key])
            else:
                state, _, _ = \
                    advance_timestep(state, params, precomputed, dt, save_tvar_history=False,
                                     total_mixing_exposure=daily_mixing_exposure)

        for field in fields(state):
            if field.name == "init_vals":
                continue
            state_history_dict[str(field.name)].append(getattr(state, field.name).clone())

    return state_history_dict, tvar_history_dict

torch_simulate_hospital_admits(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, precomputed: FluPrecomputedTensors, schedules: FluFullMetapopScheduleTensors, num_days: int, timesteps_per_day: int) -> torch.Tensor

Analogous to torch_simulate_full_history but only saves and returns hospital admits for calibration use.

Returns:

Type Description
torch.Tensor of size (num_days, L, A, R)

Returns hospital admits (the ISH to HR and HD transition variable values) for day, location, age, risk, in tensor format.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def torch_simulate_hospital_admits(state: FluFullMetapopStateTensors,
                                     params: FluFullMetapopParamsTensors,
                                     precomputed: FluPrecomputedTensors,
                                     schedules: FluFullMetapopScheduleTensors,
                                     num_days: int,
                                     timesteps_per_day: int) -> torch.Tensor:
    """
    Analogous to `torch_simulate_full_history` but only saves and
    returns hospital admits for calibration use.

    Returns:
        (torch.Tensor of size (num_days, L, A, R)):
            Returns hospital admits (the ISH to HR and HD 
            transition variable values) for day, location,
            age, risk, in tensor format.
    """

    hospital_admits_history = []

    dt = 1 / float(timesteps_per_day)

    for day in range(num_days):
        state = update_state_with_schedules(state, params, schedules, day)
        # Compute mixing exposure once per day (matching numpy metapop model)
        daily_mixing_exposure = compute_total_mixing_exposure(state, params, precomputed)
        daily_admits = None
        for timestep in range(timesteps_per_day):
            state, calibration_targets, _ = \
                advance_timestep(state, params, precomputed, dt, save_calibration_targets=True,
                                 total_mixing_exposure=daily_mixing_exposure)
            if daily_admits is None:
                daily_admits = calibration_targets["ISH_to_H"].clone()
            else:
                daily_admits = daily_admits + calibration_targets["ISH_to_H"]
        hospital_admits_history.append(daily_admits)

    return torch.stack(hospital_admits_history)

update_state_with_schedules(state: FluFullMetapopStateTensors, params: FluFullMetapopParamsTensors, schedules: FluFullMetapopScheduleTensors, day_counter: int) -> FluFullMetapopStateTensors

Returns new dataclass formed by copying the current state and updating specific values according to schedules and the simulation's current day_counter.

Returns:

Type Description
FluFullMetapopStateTensors

New state with updated schedule-related values: - flu_contact_matrix - absolute_humidity - daily_vaccines All other fields remain unchanged from the input state.

Source code in CLT_BaseModel/flu_core/flu_torch_det_components.py
def update_state_with_schedules(state: FluFullMetapopStateTensors,
                                params: FluFullMetapopParamsTensors,
                                schedules: FluFullMetapopScheduleTensors,
                                day_counter: int) -> FluFullMetapopStateTensors:
    """
    Returns new dataclass formed by copying the current `state`
    and updating specific values according to `schedules` and
    the simulation's current `day_counter`.

    Returns:
        (FluFullMetapopStateTensors):
            New state with updated schedule-related values:
              - `flu_contact_matrix`
              - `absolute_humidity`
              - `daily_vaccines`
            All other fields remain unchanged from the input `state`.
    """

    flu_contact_matrix = compute_flu_contact_matrix(params, schedules, day_counter)
    absolute_humidity = schedules.absolute_humidity[day_counter]
    daily_vaccines = schedules.daily_vaccines[day_counter]
    mobility_modifier = schedules.mobility_modifier[day_counter]

    check_and_apply_MV_reset(state, params, day_counter)

    state_new = FluFullMetapopStateTensors(
        S=state.S,
        E=state.E,
        IP=state.IP,
        ISR=state.ISR,
        ISH=state.ISH,
        IA=state.IA,
        HR=state.HR,
        HD=state.HD,
        R=state.R,
        D=state.D,
        M=state.M,
        MV=state.MV,
        absolute_humidity=absolute_humidity,
        daily_vaccines=daily_vaccines,
        flu_contact_matrix=flu_contact_matrix,
        mobility_modifier=mobility_modifier
    )

    return state_new

vaccine_preventable_events(baseline_model: MetapopModel, counterfactual_model: MetapopModel, metric_fn: Callable, **kwargs) -> float

Vaccine-preventable events = metric(baseline) βˆ’ metric(counterfactual).

A positive result means the higher-coverage counterfactual prevented events relative to the baseline.

Parameters

baseline_model : MetapopModel Model run under baseline vaccine coverage. counterfactual_model : MetapopModel Model run under a higher vaccine coverage scenario. metric_fn : callable One of cumulative_hospitalizations, cumulative_deaths, attack_rate, etc. **kwargs Forwarded to metric_fn (e.g. subpop_name, age_group).

Returns

float

Source code in CLT_BaseModel/flu_core/flu_outcomes.py
def vaccine_preventable_events(
    baseline_model: MetapopModel,
    counterfactual_model: MetapopModel,
    metric_fn: Callable,
    **kwargs,
) -> float:
    """
    Vaccine-preventable events = metric(baseline) βˆ’ metric(counterfactual).

    A positive result means the higher-coverage counterfactual prevented events
    relative to the baseline.

    Parameters
    ----------
    baseline_model : MetapopModel
        Model run under baseline vaccine coverage.
    counterfactual_model : MetapopModel
        Model run under a higher vaccine coverage scenario.
    metric_fn : callable
        One of ``cumulative_hospitalizations``, ``cumulative_deaths``,
        ``attack_rate``, etc.
    **kwargs
        Forwarded to ``metric_fn`` (e.g. ``subpop_name``, ``age_group``).

    Returns
    -------
    float
    """
    return metric_fn(baseline_model, **kwargs) - metric_fn(counterfactual_model, **kwargs)