Skip to content

Chikungunya SIR Model

A deterministic SIR model of chikungunya transmission with optional NPI and vaccination interventions.

ChikungunyaSirModel(config)

Bases: Model

Deterministic SIR model for chikungunya transmission.

Model assumptions: - Closed population with no births, deaths, or imported infections. - Homogeneous mixing within the population. - Deterministic transmission dynamics. - The epidemiological state structure is S -> I -> R. - NPI and vaccination effects are represented through modifications of the effective transmission rate rather than through additional epidemiological compartments.

The model supports two optional interventions:

  1. Non-pharmaceutical intervention (NPI)
  2. NPI effects are represented through a gradual reduction in beta.
  3. The intervention begins on simulation day 17 and reaches its maximum effect on day 28.
  4. The default maximum reduction in beta is 80%.
  5. This 80% reduction is treated as the strongest representation of the NPI intensity implemented during the outbreak.
  6. Additional NPI scenario or sensitivity analyses should therefore examine weaker intervention intensities, corresponding to reductions in beta smaller than 80%.

  7. Vaccination

  8. Vaccination eligibility is restricted to individuals aged 12 years or older.
  9. Target coverage: 40%.
  10. Vaccine efficacy against infection: 40%.
  11. Protection delay: 14 days.
  12. Vaccination starts on simulation day 0.
  13. Age eligibility is represented through the proportion of the population aged under 12 years; the S, I, and R compartments themselves are not age-stratified.

The vaccination effect follows the original R/odin model logic. Vaccination does not introduce an additional epidemiological compartment. Instead, the effective vaccination coverage is used to reduce the effective transmission rate over time.

The implementation is intended for intervention-scenario exploration in the Pandemic Simulator and preserves the mathematical structure of the original R/odin implementation.

Initialise the model and load model-specific intervention parameters.

The base Model initialiser first handles the standard framework setup, including transmission parameters, population state, and configured interventions. Model-specific NPI and vaccination parameters are then read from the Disease block of the configuration.

Source code in compartment/models/chikungunya_sir_model/model.py
def __init__(self, config):
    """
    Initialise the model and load model-specific intervention parameters.

    The base Model initialiser first handles the standard framework setup,
    including transmission parameters, population state, and configured
    interventions. Model-specific NPI and vaccination parameters are then
    read from the Disease block of the configuration.
    """

    super().__init__(config)

    disease_config = self.payload["Disease"]

    # ---------------------------------------------------------------------
    # NPI parameters
    # ---------------------------------------------------------------------
    self.npi_start_day = float(
        disease_config["npi_start_day"]
    )

    self.npi_end_day = float(
        disease_config["npi_end_day"]
    )

    self.npi_max_reduction = float(
        disease_config["npi_max_reduction"]
    )

    # ---------------------------------------------------------------------
    # Vaccination parameters
    # ---------------------------------------------------------------------
    self.prop_under12 = float(
        disease_config["prop_under12"]
    )

    self.vaccine_coverage = float(
        disease_config["vaccine_coverage"]
    )

    self.vaccine_efficacy = float(
        disease_config["vaccine_efficacy"]
    )

    self.vacc_delay = float(
        disease_config["vacc_delay"]
    )

    self.vacc_start_day = float(
        disease_config["vacc_start_day"]
    )

    self.n_days = float(
        disease_config["n_days"]
    )

prepare_initial_state()

Return the initial population matrix created by the base Model class.

The example configuration specifies the total population and the initial infected proportion. No custom age-stratified initial state is required because age eligibility is used only for vaccination coverage calculations and does not create separate age-specific compartments.

Source code in compartment/models/chikungunya_sir_model/model.py
def prepare_initial_state(self):
    """
    Return the initial population matrix created by the base Model class.

    The example configuration specifies the total population and the
    initial infected proportion. No custom age-stratified initial state is
    required because age eligibility is used only for vaccination coverage
    calculations and does not create separate age-specific compartments.
    """

    return self.population_matrix

custom_npi(t, beta)

Apply a gradual NPI-related reduction in the transmission rate beta.

The intervention begins on npi_start_day and increases linearly until npi_end_day, after which the maximum reduction is maintained.

In the default configuration: - npi_start_day = 17 - npi_end_day = 28 - npi_max_reduction = 0.80

The 80% reduction is treated as the strongest representation of the NPI intensity implemented during the outbreak. Therefore, additional NPI scenario or sensitivity analyses should represent weaker interventions by specifying reductions in beta below 80%.

If the NPI intervention is disabled, the baseline beta is returned unchanged.

Source code in compartment/models/chikungunya_sir_model/model.py
def custom_npi(self, t, beta):
    """
    Apply a gradual NPI-related reduction in the transmission rate beta.

    The intervention begins on npi_start_day and increases linearly
    until npi_end_day, after which the maximum reduction is maintained.

    In the default configuration:
    - npi_start_day = 17
    - npi_end_day = 28
    - npi_max_reduction = 0.80

    The 80% reduction is treated as the strongest representation of
    the NPI intensity implemented during the outbreak. Therefore,
    additional NPI scenario or sensitivity analyses should represent
    weaker interventions by specifying reductions in beta below 80%.

    If the NPI intervention is disabled, the baseline beta is returned
    unchanged.
    """

    # NPI OFF:
    # return the baseline transmission rate unchanged.
    if not self._intervention_enabled("npi"):
        return beta

    # Duration of the linear NPI ramp.
    duration = (
        self.npi_end_day
        - self.npi_start_day
    )

    # Avoid division by zero if start and end days were ever configured
    # to be identical.
    safe_duration = jnp.maximum(
        duration,
        1e-10,
    )

    # Linear progress of the NPI effect:
    #
    # t < start day     -> progress = 0
    # start <= t < end  -> progress increases from 0 to 1
    # t >= end day      -> progress = 1
    progress = (
        (t - self.npi_start_day)
        / safe_duration
    )

    progress = jnp.clip(
        progress,
        0.0,
        1.0,
    )

    # Current proportional reduction in beta.
    reduction = (
        self.npi_max_reduction
        * progress
    )

    # Effective transmission rate after applying the NPI.
    effective_beta = (
        beta
        * (1.0 - reduction)
    )

    return effective_beta

custom_vaccination(t, beta)

Apply vaccination following the original R/odin model logic.

Vaccination is represented through a time-varying reduction in beta. No separate vaccinated epidemiological compartment is introduced.

Age eligibility is handled using prop_under12: only the population aged 12 years or older is considered eligible, but S, I, and R are not explicitly age-stratified.

Original R/odin logic:

effective_vacc_days = n_days - vacc_start_day - vacc_delay

can_vaccinate = effective_vacc_days > 0

daily_vacc_prop = if can_vaccinate: vaccine_coverage / effective_vacc_days else: 0

effective_days = max( 0, t - (vacc_start_day + vacc_delay) + 1 )

vaccine_cov_current = min( vaccine_coverage, daily_vacc_prop * effective_days * vaccine_efficacy )

vaccine_cov_total_population = vaccine_cov_current * (1 - prop_under12)

transmission_reduction = 1 - vaccine_cov_total_population

effective_beta = beta * transmission_reduction

The placement of vaccine efficacy inside vaccine_cov_current and the +1 term in effective_days are retained exactly from the original R/odin implementation to preserve the original model dynamics.

Source code in compartment/models/chikungunya_sir_model/model.py
def custom_vaccination(self, t, beta):
    """
    Apply vaccination following the original R/odin model logic.

    Vaccination is represented through a time-varying reduction in beta.
    No separate vaccinated epidemiological compartment is introduced.

    Age eligibility is handled using prop_under12:
    only the population aged 12 years or older is considered eligible,
    but S, I, and R are not explicitly age-stratified.

    Original R/odin logic:

    effective_vacc_days =
        n_days - vacc_start_day - vacc_delay

    can_vaccinate =
        effective_vacc_days > 0

    daily_vacc_prop =
        if can_vaccinate:
            vaccine_coverage / effective_vacc_days
        else:
            0

    effective_days =
        max(
            0,
            t - (vacc_start_day + vacc_delay) + 1
        )

    vaccine_cov_current =
        min(
            vaccine_coverage,
            daily_vacc_prop *
            effective_days *
            vaccine_efficacy
        )

    vaccine_cov_total_population =
        vaccine_cov_current *
        (1 - prop_under12)

    transmission_reduction =
        1 - vaccine_cov_total_population

    effective_beta =
        beta * transmission_reduction

    The placement of vaccine efficacy inside vaccine_cov_current and the
    +1 term in effective_days are retained exactly from the original
    R/odin implementation to preserve the original model dynamics.
    """

    # ---------------------------------------------------------------------
    # Vaccination OFF
    # ---------------------------------------------------------------------
    # When vaccination is disabled, this is equivalent to setting
    # vaccination coverage to zero in the original R analysis.
    #
    # beta is therefore returned unchanged.
    if not self._intervention_enabled("vaccination"):
        return beta

    # ---------------------------------------------------------------------
    # Number of days available for effective vaccination
    # ---------------------------------------------------------------------
    # Original R:
    #
    # effective_vacc_days <-
    #     n_days - vacc_start_day - vacc_delay
    #
    # The protection delay is subtracted because vaccine effects are not
    # assumed to begin immediately after vaccination starts.
    effective_vacc_days = (
        self.n_days
        - self.vacc_start_day
        - self.vacc_delay
    )

    can_vaccinate = (
        effective_vacc_days > 0.0
    )

    # Numerical safeguard for the denominator.
    #
    # This does not alter the intended logic because daily_vacc_prop is
    # set to zero through jnp.where when effective_vacc_days <= 0.
    safe_effective_vacc_days = jnp.maximum(
        effective_vacc_days,
        1e-10,
    )

    # ---------------------------------------------------------------------
    # Daily vaccination proportion
    # ---------------------------------------------------------------------
    # Original R:
    #
    # daily_vacc_prop <-
    #     if (can_vaccinate)
    #         vacc_coverage_adult_target /
    #         effective_vacc_days
    #     else 0
    #
    # This spreads the target vaccination coverage uniformly across the
    # available vaccination period.
    daily_vacc_prop = jnp.where(
        can_vaccinate,
        self.vaccine_coverage
        / safe_effective_vacc_days,
        0.0,
    )

    # ---------------------------------------------------------------------
    # Number of days for which vaccine protection has become effective
    # ---------------------------------------------------------------------
    # Original R:
    #
    # effective_days <-
    #     max(
    #         0,
    #         t - (vacc_start_day + vacc_delay) + 1
    #     )
    #
    # The +1 follows the original R/odin implementation exactly and is
    # retained to preserve its timing convention.
    effective_days = jnp.maximum(
        0.0,
        t
        - (
            self.vacc_start_day
            + self.vacc_delay
        )
        + 1.0,
    )

    # ---------------------------------------------------------------------
    # Current effective vaccination coverage among the eligible population
    # ---------------------------------------------------------------------
    # Original R:
    #
    # vaccine_cov_current <-
    #     min(
    #         vacc_coverage_adult_target,
    #         daily_vacc_prop *
    #         effective_days *
    #         vaccine_efficacy
    #     )
    #
    # Vaccine efficacy is intentionally multiplied within this effective
    # coverage term, exactly as in the original R/odin implementation.
    #
    # The current effective coverage cannot exceed the configured target
    # vaccination coverage.
    vaccine_cov_current = jnp.minimum(
        self.vaccine_coverage,
        daily_vacc_prop
        * effective_days
        * self.vaccine_efficacy,
    )

    # ---------------------------------------------------------------------
    # Convert effective coverage among age 12+ to the total population
    # ---------------------------------------------------------------------
    # Original R:
    #
    # vaccine_cov_total_population <-
    #     vaccine_cov_current *
    #     (1 - prop_under12)
    #
    # This is how the age restriction is represented without explicitly
    # creating age-stratified S, I, and R compartments.
    vaccine_cov_total_population = (
        vaccine_cov_current
        * (1.0 - self.prop_under12)
    )

    # ---------------------------------------------------------------------
    # Remaining transmission fraction after vaccination
    # ---------------------------------------------------------------------
    # Original R:
    #
    # transmission_reduction <-
    #     1 - vaccine_cov_total_population
    #
    # Despite the original variable name, this quantity is the remaining
    # fraction of baseline transmission after accounting for vaccination.
    transmission_reduction = (
        1.0
        - vaccine_cov_total_population
    )

    # ---------------------------------------------------------------------
    # Effective transmission rate after vaccination
    # ---------------------------------------------------------------------
    # Original R:
    #
    # effective_beta <-
    #     sir_beta * transmission_reduction
    #
    # Vaccination affects infection dynamics only through effective_beta.
    effective_beta = (
        beta
        * transmission_reduction
    )

    return effective_beta

equation(y, t, p)

Compute the compartment derivatives for one integration step.

The derivatives are written explicitly rather than relying entirely on the framework's automatic equation builder. This preserves the custom force of infection in which beta is modified sequentially by the optional NPI and vaccination mechanisms.

Parameters:

Name Type Description Default
y

Current compartment values in framework-defined order.

required
t

Current simulation time in days since the simulation start.

required
p

Packed model parameter vector.

required

Returns:

Type Description

JAX array containing derivatives in self.compartment_list order.

Source code in compartment/models/chikungunya_sir_model/model.py
def equation(self, y, t, p):
    """
    Compute the compartment derivatives for one integration step.

    The derivatives are written explicitly rather than relying entirely
    on the framework's automatic equation builder. This preserves the
    custom force of infection in which beta is modified sequentially by
    the optional NPI and vaccination mechanisms.

    Args:
        y:
            Current compartment values in framework-defined order.
        t:
            Current simulation time in days since the simulation start.
        p:
            Packed model parameter vector.

    Returns:
        JAX array containing derivatives in self.compartment_list order.
    """

    C = self.COMPARTMENTS

    # ---------------------------------------------------------------------
    # Parameters
    # ---------------------------------------------------------------------
    # Convert the framework parameter vector into a named dictionary.
    params = self._unpack_params(p)

    # ---------------------------------------------------------------------
    # Current states
    # ---------------------------------------------------------------------
    # Never hardcode compartment positions. Mapping through
    # self.compartment_list preserves compatibility with the framework's
    # automatically generated cumulative compartments.
    states = {
        c: y[i]
        for i, c in enumerate(self.compartment_list)
    }

    S = states[C.S]
    I = states[C.I]
    R = states[C.R]

    # ---------------------------------------------------------------------
    # Population total
    # ---------------------------------------------------------------------
    # In the original R model, population_total is fixed.
    #
    # Here:
    #
    # dS/dt + dI/dt + dR/dt = 0
    #
    # so S + I + R remains constant and is equal to the initial total
    # population throughout the simulation.
    #
    # This preserves the fixed population_total assumption of the
    # original R/odin implementation.
    population_total = (
        S
        + I
        + R
    )

    # ---------------------------------------------------------------------
    # Baseline beta
    # ---------------------------------------------------------------------
    beta = params["beta"]

    # ---------------------------------------------------------------------
    # Optional NPI
    # ---------------------------------------------------------------------
    # The custom NPI modifies beta only when the NPI intervention is
    # enabled.
    #
    # In the default NPI scenario, beta is gradually reduced by up to
    # 80%, representing the strongest implemented NPI intensity
    # considered in this model.
    #
    # Weaker NPI scenarios should use smaller reductions in beta.
    #
    # For the vaccine-only example configuration, NPI is not enabled,
    # so custom_npi() simply returns the baseline beta.
    beta = self.custom_npi(
        t,
        beta,
    )

    # ---------------------------------------------------------------------
    # Optional vaccination
    # ---------------------------------------------------------------------
    # Vaccination is applied after the optional NPI modification.
    #
    # In a vaccine-only run, beta entering this function is the baseline
    # beta. In a combined NPI + vaccination run, beta already contains
    # the NPI-related reduction.
    effective_beta = self.custom_vaccination(
        t,
        beta,
    )

    # gamma has already been converted by the framework from the
    # configured DAYS representation to the corresponding recovery rate.
    gamma = params["gamma"]

    # ---------------------------------------------------------------------
    # Force of infection
    # ---------------------------------------------------------------------
    # Original R:
    #
    # s_to_i <-
    #     effective_beta *
    #     state_s *
    #     state_i /
    #     population_total
    #
    # This is a frequency-dependent SIR force of infection.
    s_to_i = (
        effective_beta
        * S
        * I
        / (
            population_total
            + 1e-10
        )
    )

    # ---------------------------------------------------------------------
    # Recovery flow
    # ---------------------------------------------------------------------
    i_to_r = (
        gamma
        * I
    )

    # ---------------------------------------------------------------------
    # Initialise derivatives
    # ---------------------------------------------------------------------
    # self.compartment_list also includes framework-generated cumulative
    # compartments such as I_total and R_total.
    derivs = {
        c: jnp.zeros_like(I)
        for c in self.compartment_list
    }

    # ---------------------------------------------------------------------
    # Epidemiological SIR dynamics
    # ---------------------------------------------------------------------
    derivs[C.S] = (
        -s_to_i
    )

    derivs[C.I] = (
        s_to_i
        - i_to_r
    )

    derivs[C.R] = (
        i_to_r
    )

    # ---------------------------------------------------------------------
    # Framework cumulative totals
    # ---------------------------------------------------------------------
    # I_total records cumulative incident infections.
    # R_total records cumulative recoveries.
    #
    # These are bookkeeping/output states only and do not contribute
    # to the force of infection or active population total.
    if f"{C.I}_total" in derivs:
        derivs[
            f"{C.I}_total"
        ] = s_to_i

    if f"{C.R}_total" in derivs:
        derivs[
            f"{C.R}_total"
        ] = i_to_r

    # ---------------------------------------------------------------------
    # Preserve framework-defined compartment order
    # ---------------------------------------------------------------------
    # The framework expects derivatives to be returned in exactly the
    # same order as self.compartment_list.
    return jnp.stack(
        [
            derivs[c]
            for c in self.compartment_list
        ]
    )