26623145945

Committed 29 May 2026 07:00AM UTC coverage: 90.685% (-0.2%) from 90.92%

Build # 26623145945

Build Type

push

github

Committed by

web-flow

Commit Message

ci: migrate from Cirrus to GitHub Actions (#76)

* chore(ci): register slow marker and add py311/py312 tox envs

Pre-flight hygiene for the CI migration off Cirrus. Two provider-independent
fixes that de-risk the upcoming matrix:

- Register the slow pytest marker in setup.cfg so CI logs no longer carry a
  PytestUnknownMarkWarning on every run. The slow tests already exist
  (KirmansAnts stationary KS, two common tests); only the registration was
  missing.
- Add explicit py311 and py312 tox envs mirroring the existing py313 env.
  Without these, tox -e py311 falls through to the base [testenv] and runs
  the full slow suite, which would surprise the matrix CI.

* ci: migrate from Cirrus to GitHub Actions

Cirrus-CI shuts down 2026-05-31. Replace it with a 2-OS x 3-Python matrix
on GitHub Actions, restoring the wheel-install pattern and rewiring the
README badges. No-op for local development.

Workflow shape:
- build job (ubuntu-latest, Python 3.13) runs 'tox -e clean,build' and
  uploads dist/* as an artifact, so every matrix leg installs the SAME
  wheel and catches packaging bugs that per-job rebuilds miss.
- test matrix (ubuntu-latest, windows-latest) x (3.11, 3.12, 3.13) with
  fail-fast: false to mirror Cirrus' independent-task semantics.
- actions/checkout uses fetch-depth: 0 because setuptools_scm is
  configured with version_scheme = 'no-guess-dev' and refuses to build
  on shallow clones without tag history.
- actions/setup-python@v5 caches pip keyed on pyproject.toml + uv.lock.
  The torch CPU wheel is ~200 MB; the cache is the single biggest perf
  win across the matrix.
- Windows enables long-path support via PowerShell before any install,
  porting the Cirrus registry tweak.
- Coverage uploads per-leg with coverallsapp/github-action@v2
  parallel: true; a final coverage-finish job calls parallel-finished
  so Coveralls merges the six legs into one report.
- schedule cron at 04:00 UTC nightly runs the FULL suite (slow + fast)
  on master... (continued)

Coverage Stats

495 of 538 branches covered (92.01%)

Branch coverage included in aggregate %.

1 of 1 new or added line in 1 file covered. (100.0%)

7 existing lines in 2 files now uncovered.

3331 of 3681 relevant lines covered (90.49%)

5.43 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

93.21

/src/macrostat/core/behavior.py

"""
Behavior classes for the MacroStat model.
"""

__author__ = ["Karl Naumann-Woleske"]
__credits__ = ["Karl Naumann-Woleske"]
__license__ = "MIT"
__maintainer__ = ["Karl Naumann-Woleske"]

import logging

import numpy as np
import torch

from macrostat.core.parameters import Parameters
from macrostat.core.scenarios import Scenarios
from macrostat.core.variables import Variables

logger = logging.getLogger(__name__)


class Behavior(torch.nn.Module):
    """Base class for the behavior of the MacroStat model."""

    supports_differentiable: bool = True

    def __init__(
        self,
        parameters: Parameters,
        scenarios: Scenarios,
        variables: Variables,
        scenario: int = 0,
        differentiable: bool = False,
        debug: bool = False,
    ):
        """Initialize the behavior of the MacroStat model.

        Parameters
        ----------
        parameters: macrostat.core.parameters.Parameters
            The parameters of the model.
        scenarios: macrostat.core.scenarios.Scenarios
            The scenarios of the model.
        variables: macrostat.core.variables.Variables
            The variables of the model.
        scenario: int
            The scenario to use for the model run.
        debug: bool
            Whether to print debug information.

        Raises
        ------
        RuntimeError
            If ``differentiable=True`` is requested for a subclass that sets
            ``supports_differentiable = False`` (e.g. SDE models with
            non-differentiable rejection sampling).
        """
        # Initialize the parent class
        super().__init__()

        if differentiable and not self.supports_differentiable:
            raise RuntimeError(
                f"{type(self).__name__} sets supports_differentiable=False "
                "and cannot be constructed with differentiable=True."
            )

        # Initialize the parameters. Keep the Parameters instance itself
        # so step-time calls can always ask for a fresh resolver and
        # constraint list. Snapshotting them at init time would freeze
        # the view of the parameter layout at construction and make any
        # future mutation API silently desynchronise.
        self.parameters = parameters
        self.params = parameters.to_nn_parameters()
        self.hyper = parameters.hyper

        # Initialize the scenarios
        self.scenarios = scenarios.to_nn_parameters(scenario=scenario)
        self.scenarioID = scenario

        # Initialize the variables
        self.variables = variables

        # Settings
        self.differentiable = differentiable
        self.debug = debug

        # Per-instance numpy RNG seeded from hyperparameter "seed". Used by
        # subclasses that need a numpy.random.Generator without touching the
        # global numpy state — required for parallel runs (e.g. Hessian
        # computation across parameter perturbations).
        self.numpy_rng: np.random.Generator = np.random.default_rng(self.hyper["seed"])

    ############################################################################
    # Simulation of the model
    ############################################################################

    def forward(self):
        """Forward pass of the behavior.

        This should include the model's main loop, and is implemented as a placeholder.
        The idea is for users to implement an initialize() and step() function,
        which will be called by the forward() function.

        If there are additional steps necessary, users may wish to overwrite this function.
        """
        # Set the seed
        torch.manual_seed(self.hyper["seed"])
        self.numpy_rng = np.random.default_rng(self.hyper["seed"])

        # Initialize the output tensors
        self.state, self.history = self.variables.initialize_tensors()

        # Initialize the model
        logger.debug(
            f"Initializing model (t=0...{self.hyper['timesteps_initialization']})"
        )
        self.initialize()

        for t in range(self.hyper["timesteps_initialization"] + 1):
            self.variables.record_state(t, self.state)

        for t in range(self.hyper["timesteps_initialization"] + 1):
            self.history = self.variables.update_history(self.state)

        # Initialize the prior and state
        self.prior = self.state

        # Run the model for the remaining timesteps
        logger.debug(
            f"Simulating model (t={self.hyper['timesteps_initialization'] + 1}...{self.hyper['timesteps']})"
        )

        for t in range(self.hyper["timesteps_initialization"], self.hyper["timesteps"]):
            self.state = self.variables.new_state()
            # Get scenario series for this point in time
            idx = torch.where(
                torch.arange(self.hyper["timesteps"]) == t,
                torch.ones(1),
                torch.zeros(1),
            )
            scenario = {k: idx @ v for k, v in self.scenarios.items()}

            # Apply parameter shocks
            params = self.apply_parameter_shocks(t, scenario)

            # Step the model
            self.step(t=t, scenario=scenario, params=params)

            # Store the outputs
            self.variables.record_state(t, self.state)
            self.history = self.variables.update_history(self.state)
            self.prior = self.state

        return self.variables.gather_timeseries()

    def initialize(self):
        """Initialize the behavior.

        This should include the model's initialization steps, and set all of the
        necessary state variables. They only need to be set for one period, and
        will then be copied to the history and prior to be used in the step function.
        """
        raise NotImplementedError("Behavior.initialize() to be implemented by model")

    def step(self, t: int, scenario: dict, params: dict | None = None):
        """Step function of the behavior.

        This should include the model's main loop.

        Parameters
        ----------
        t: int
            The current timestep.
        scenario: dict
            The scenario information for the current timestep.
        """
        raise NotImplementedError("Behavior.step() to be implemented by model")

    def apply_parameter_shocks(self, t: int, scenario: dict):
        """Apply parameter shocks to the model.

        Any parameter in the model can be shocked/changed during the simulation
        using the scenario information. Specifically, for a parameter alpha, the
        user can pass two types of potential shocks:
        1. An multiplicative shock, generically named alpha_multiply
        2. An additive shock, generically named alpha_add

        This function will apply the shocks to the parameters, and return a
        dictionary with the updated parameters. The application of the shocks is
        independent, that is, the multiplicative shock does not affect the additive
        shock and vice versa. This is done by first applying the multiplicative
        shock, and then the additive shock.

        Parameters
        ----------
        t: int
            The current timestep.
        scenario: dict
            The scenario information for the current timestep.

        Returns
        -------
        dict
            A dictionary with the updated parameters.
        """
        # Optional sectoral structure for vector/matrix parameters
        vsecs = self.hyper.get("vector_sectors", [])
        n = len(vsecs)
        if n > 0:
            one, zero = torch.ones(n), torch.zeros(n)
            sec_vectors = {
                s: torch.where(torch.arange(n) == i, one, zero)
                for i, s in enumerate(vsecs)
            }

            # Generate index matrices per sector pair
            pairs = [(row, col) for row in vsecs for col in vsecs]
            one, zero = torch.ones(n, n), torch.zeros(n, n)
            sec_matrices = {
                s: torch.where(torch.arange(n * n).reshape(n, n) == i, one, zero)
                for i, s in enumerate(pairs)
            }
        else:
            sec_vectors = {}
            sec_matrices = {}

        params = {}
        for key, value in self.params.items():
            if len(value.shape) == 0:
                mul, add = torch.tensor(1.0), torch.tensor(0.0)
                if f"{key}_multiply" in scenario:
                    mul = scenario[f"{key}_multiply"]
                if f"{key}_add" in scenario:
                    add = scenario[f"{key}_add"]

            else:
                add = torch.zeros_like(value)
                mul = torch.ones_like(value)

                if len(value.shape) == 1 and sec_vectors:
                    for s, ix in sec_vectors.items():
                        if f"{s}_{key}_multiply" in scenario:
                            mul = mul * (ix * scenario[f"{s}_{key}_multiply"])
                        if f"{s}_{key}_add" in scenario:
                            add = add + (ix * scenario[f"{s}_{key}_add"])

                elif len(value.shape) == 2 and sec_matrices:
                    for rowcol, ix in sec_matrices.items():
                        s = f"{rowcol[0]}_{rowcol[1]}"

                        if f"{s}_{key}_multiply" in scenario:
                            mul = mul * (ix * scenario[f"{s}_{key}_multiply"])
                        if f"{s}_{key}_add" in scenario:
                            add = add + (ix * scenario[f"{s}_{key}_add"])

            params[key] = value * mul + add

        # Enforce adding-up constraints (differentiable). The resolver
        # and constraint list are fetched fresh on every call so the
        # step-time view of the parameter layout always matches the
        # current Parameters instance.
        constraints = self.parameters.get_constraints()
        if constraints:
            resolver = self.parameters.get_constraint_resolver()
            for c in constraints:
                c.apply(params, resolver)

        return params

    ############################################################################
    # Steady State
    ############################################################################

    def compute_theoretical_steady_state(self, **kwargs):
        """Compute the theoretical steady state of the model.

        This process generally follows the structure of the forward() function,
        but instead of simulating the model, the steady state is computed at
        each timestep. Therefore, (1) the model is initialized, and (2) for
        each timestep the parameter and scenario information is passed to the
        compute_theoretical_steady_state_per_step() function that computes the
        steady state at that timestep.

        Parameters
        ----------
        **kwargs: dict
            Additional keyword arguments.
        """
        # Set the seed
        torch.manual_seed(self.hyper["seed"])
        self.numpy_rng = np.random.default_rng(self.hyper["seed"])

        # Initialize the output tensors
        self.state, self.history = self.variables.initialize_tensors()

        # Initialize the model
        info = f"(t=0...{self.hyper['timesteps_initialization']})"
        logger.debug(f"Initializing model {info}")
        self.initialize()

        for t in range(self.hyper["timesteps_initialization"] + 1):
            self.variables.record_state(t, self.state)

        for t in range(self.hyper["timesteps_initialization"] + 1):
            self.history = self.variables.update_history(self.state)

        # Initialize prior from the initialization state
        self.prior = self.state

        # Compute the steady state
        info = f"(t={self.hyper['timesteps_initialization'] + 1}...{self.hyper['timesteps']})"
        logger.debug(f"Computing theoretical steady state {info}")

        for t in range(
            self.hyper["timesteps_initialization"] + 1, self.hyper["timesteps"]
        ):
            self.state = self.variables.new_state()

            # Get scenario series for this point in time
            idx = torch.where(
                torch.arange(self.hyper["timesteps"]) == t,
                torch.ones(1),
                torch.zeros(1),
            )
            scenario = {k: idx @ v for k, v in self.scenarios.items()}

            # Apply parameter shocks
            params = self.apply_parameter_shocks(t, scenario)

            # Compute the steady state
            self.compute_theoretical_steady_state_per_step(
                t=t, params=params, scenario=scenario
            )

            # Store the outputs
            self.variables.record_state(t, self.state)
            self.history = self.variables.update_history(self.state)
            self.prior = self.state

        return None

    def compute_theoretical_steady_state_per_step(self, **kwargs):
        """Compute the theoretical steady state of the model per step."""
        raise NotImplementedError(
            "Behavior.compute_theoretical_steady_state_per_step() to be implemented by model"
        )

    ############################################################################
    # Some Differentiable PyTorch Alternatives
    ############################################################################

    def diffwhere(self, condition, x1, x2):
        """Where condition that is differentiable with respect to the condition.

        Requires:
            self.hyper['diffwhere'] = True
            self.hyper['sigmoid_constant'] as a large number

        Note: For non-NaN/inf, where(x > eps, z, y) is (x - eps > 0) * (z - y) + y,
        so we can use the sigmoid function to approximate the where function.

        Parameters
        ----------
        condition : torch.Tensor
            Condition to be evaluated expressed as x - eps
        x1 : torch.Tensor
            Value to be returned if condition is True
        x2 : torch.Tensor
            Value to be returned if condition is False
        """
        sig = torch.sigmoid(torch.mul(condition, self.hyper["sigmoid_constant"]))
        return torch.add(torch.mul(sig, torch.sub(x1, x2)), x2)

    def tanhmask(self, x):
        """Convert a variable into 0 (x<0) and 1 (x>0)

        Requires:
            self.hyper['tanh_constant'] as a large number

        Parameters
        ----------
        x: torch.Tensor
            The variable to be converted.

        """
        kwg = {"dtype": torch.float64, "requires_grad": True}
        return torch.div(
            torch.add(
                torch.ones(x.size(), **kwg),
                torch.tanh(torch.mul(x, self.hyper["tanh_constant"])),
            ),
            torch.tensor(2.0, **kwg),
        )

    def diffmin(self, x1, x2):
        """Smooth approximation to the minimum
        B: https://mathoverflow.net/questions/35191/a-differentiable-approximation-to-the-minimum-function

        Requires:
            self.hyper['min_constant'] as a large number

        Parameters
        ----------
        x1: torch.Tensor
            The first variable to be compared.
        x2: torch.Tensor
            The second variable to be compared.
        """
        r = self.hyper["min_constant"]
        pt1 = torch.exp(torch.mul(x1, -1 * r))
        pt2 = torch.exp(torch.mul(x2, -1 * r))
        return torch.mul(-1 / r, torch.log(torch.add(pt1, pt2)))

    def diffmax(self, x1, x2):
        """Smooth approximation to the minimum
        B: https://mathoverflow.net/questions/35191/a-differentiable-approximation-to-the-minimum-function

        Requires:
            self.hyper['max_constant'] as a large number

        Parameters
        ----------
        x1: torch.Tensor
            The first variable to be compared.
        x2: torch.Tensor
            The second variable to be compared.
        """
        r = self.hyper["max_constant"]
        pt1 = torch.exp(torch.mul(x1, r))
        pt2 = torch.exp(torch.mul(x2, r))
        return torch.mul(1 / r, torch.log(torch.add(pt1, pt2)))

    def diffmin_v(self, x):
        """Smooth approximation to the minimum. See diffmin

        Parameters
        ----------
        x: torch.Tensor
            The variable to be converted.

        Requires:
            self.hyper['min_constant'] as a large number
        """
        r = self.hyper["min_constant"]
        temp = torch.exp(torch.mul(x, -1 * r))
        return torch.mul(-1 / r, torch.log(torch.sum(temp)))

    def diffmax_v(self, x):
        """Smooth approximation to the maximum for a tensor. See diffmax

        Requires:
            self.hyper['max_constant'] as a large number

        Parameters
        ----------
        x: torch.Tensor
            The variable to be converted.
        """
        r = self.hyper["max_constant"]
        temp = torch.exp(torch.mul(x, r))
        return torch.mul(1 / r, torch.log(torch.sum(temp)))


if __name__ == "__main__":
    pass

1	"""
2	Behavior classes for the MacroStat model.
3	"""
4
5	__author__ = ["Karl Naumann-Woleske"]	6✔
6	__credits__ = ["Karl Naumann-Woleske"]	6✔
7	__license__ = "MIT"	6✔
8	__maintainer__ = ["Karl Naumann-Woleske"]	6✔
9
10	import logging	6✔
11
12	import numpy as np	6✔
13	import torch	6✔
14
15	from macrostat.core.parameters import Parameters	6✔
16	from macrostat.core.scenarios import Scenarios	6✔
17	from macrostat.core.variables import Variables	6✔
18
19	logger = logging.getLogger(__name__)	6✔
20
21
22	class Behavior(torch.nn.Module):	6✔
23	"""Base class for the behavior of the MacroStat model."""
24
25	supports_differentiable: bool = True	6✔
26
27	def __init__(	6✔
28	self,
29	parameters: Parameters,
30	scenarios: Scenarios,
31	variables: Variables,
32	scenario: int = 0,
33	differentiable: bool = False,
34	debug: bool = False,
35	):
36	"""Initialize the behavior of the MacroStat model.
37
38	Parameters
39	----------
40	parameters: macrostat.core.parameters.Parameters
41	The parameters of the model.
42	scenarios: macrostat.core.scenarios.Scenarios
43	The scenarios of the model.
44	variables: macrostat.core.variables.Variables
45	The variables of the model.
46	scenario: int
47	The scenario to use for the model run.
48	debug: bool
49	Whether to print debug information.
50
51	Raises
52	------
53	RuntimeError
54	If ``differentiable=True`` is requested for a subclass that sets
55	``supports_differentiable = False`` (e.g. SDE models with
56	non-differentiable rejection sampling).
57	"""
58	# Initialize the parent class
59	super().__init__()	6✔
60
61	if differentiable and not self.supports_differentiable:	6✔
62	raise RuntimeError(	×
63	f"{type(self).__name__} sets supports_differentiable=False "
64	"and cannot be constructed with differentiable=True."
65	)
66
67	# Initialize the parameters. Keep the Parameters instance itself
68	# so step-time calls can always ask for a fresh resolver and
69	# constraint list. Snapshotting them at init time would freeze
70	# the view of the parameter layout at construction and make any
71	# future mutation API silently desynchronise.
72	self.parameters = parameters	6✔
73	self.params = parameters.to_nn_parameters()	6✔
74	self.hyper = parameters.hyper	6✔
75
76	# Initialize the scenarios
77	self.scenarios = scenarios.to_nn_parameters(scenario=scenario)	6✔
78	self.scenarioID = scenario	6✔
79
80	# Initialize the variables
81	self.variables = variables	6✔
82
83	# Settings
84	self.differentiable = differentiable	6✔
85	self.debug = debug	6✔
86
87	# Per-instance numpy RNG seeded from hyperparameter "seed". Used by
88	# subclasses that need a numpy.random.Generator without touching the
89	# global numpy state — required for parallel runs (e.g. Hessian
90	# computation across parameter perturbations).
91	self.numpy_rng: np.random.Generator = np.random.default_rng(self.hyper["seed"])	6✔
92
93	############################################################################
94	# Simulation of the model
95	############################################################################
96
97	def forward(self):	6✔
98	"""Forward pass of the behavior.
99
100	This should include the model's main loop, and is implemented as a placeholder.
101	The idea is for users to implement an initialize() and step() function,
102	which will be called by the forward() function.
103
104	If there are additional steps necessary, users may wish to overwrite this function.
105	"""
106	# Set the seed
107	torch.manual_seed(self.hyper["seed"])	6✔
108	self.numpy_rng = np.random.default_rng(self.hyper["seed"])	6✔
109
110	# Initialize the output tensors
111	self.state, self.history = self.variables.initialize_tensors()	6✔
112
113	# Initialize the model
114	logger.debug(	6✔
115	f"Initializing model (t=0...{self.hyper['timesteps_initialization']})"
116	)
117	self.initialize()	6✔
118
119	for t in range(self.hyper["timesteps_initialization"] + 1):	6✔
120	self.variables.record_state(t, self.state)	6✔
121
122	for t in range(self.hyper["timesteps_initialization"] + 1):	6✔
123	self.history = self.variables.update_history(self.state)	6✔
124
125	# Initialize the prior and state
126	self.prior = self.state	6✔
127
128	# Run the model for the remaining timesteps
129	logger.debug(	6✔
130	f"Simulating model (t={self.hyper['timesteps_initialization'] + 1}...{self.hyper['timesteps']})"
131	)
132
133	for t in range(self.hyper["timesteps_initialization"], self.hyper["timesteps"]):	6✔
134	self.state = self.variables.new_state()	6✔
135	# Get scenario series for this point in time
136	idx = torch.where(	6✔
137	torch.arange(self.hyper["timesteps"]) == t,
138	torch.ones(1),
139	torch.zeros(1),
140	)
141	scenario = {k: idx @ v for k, v in self.scenarios.items()}	6✔
142
143	# Apply parameter shocks
144	params = self.apply_parameter_shocks(t, scenario)	6✔
145
146	# Step the model
147	self.step(t=t, scenario=scenario, params=params)	6✔
148
149	# Store the outputs
150	self.variables.record_state(t, self.state)	6✔
151	self.history = self.variables.update_history(self.state)	6✔
152	self.prior = self.state	6✔
153
154	return self.variables.gather_timeseries()	6✔
155
156	def initialize(self):	6✔
157	"""Initialize the behavior.
158
159	This should include the model's initialization steps, and set all of the
160	necessary state variables. They only need to be set for one period, and
161	will then be copied to the history and prior to be used in the step function.
162	"""
163	raise NotImplementedError("Behavior.initialize() to be implemented by model")
164
165	def step(self, t: int, scenario: dict, params: dict \| None = None):	6✔
166	"""Step function of the behavior.
167
168	This should include the model's main loop.
169
170	Parameters
171	----------
172	t: int
173	The current timestep.
174	scenario: dict
175	The scenario information for the current timestep.
176	"""
177	raise NotImplementedError("Behavior.step() to be implemented by model")
178
179	def apply_parameter_shocks(self, t: int, scenario: dict):	6✔
180	"""Apply parameter shocks to the model.
181
182	Any parameter in the model can be shocked/changed during the simulation
183	using the scenario information. Specifically, for a parameter alpha, the
184	user can pass two types of potential shocks:
185	1. An multiplicative shock, generically named alpha_multiply
186	2. An additive shock, generically named alpha_add
187
188	This function will apply the shocks to the parameters, and return a
189	dictionary with the updated parameters. The application of the shocks is
190	independent, that is, the multiplicative shock does not affect the additive
191	shock and vice versa. This is done by first applying the multiplicative
192	shock, and then the additive shock.
193
194	Parameters
195	----------
196	t: int
197	The current timestep.
198	scenario: dict
199	The scenario information for the current timestep.
200
201	Returns
202	-------
203	dict
204	A dictionary with the updated parameters.
205	"""
206	# Optional sectoral structure for vector/matrix parameters
207	vsecs = self.hyper.get("vector_sectors", [])	6✔
208	n = len(vsecs)	6✔
209	if n > 0:	6✔
210	one, zero = torch.ones(n), torch.zeros(n)	6✔
211	sec_vectors = {	6✔
212	s: torch.where(torch.arange(n) == i, one, zero)
213	for i, s in enumerate(vsecs)
214	}
215
216	# Generate index matrices per sector pair
217	pairs = [(row, col) for row in vsecs for col in vsecs]	6✔
218	one, zero = torch.ones(n, n), torch.zeros(n, n)	6✔
219	sec_matrices = {	6✔
220	s: torch.where(torch.arange(n * n).reshape(n, n) == i, one, zero)
221	for i, s in enumerate(pairs)
222	}
223	else:
224	sec_vectors = {}	6✔
225	sec_matrices = {}	6✔
226
227	params = {}	6✔
228	for key, value in self.params.items():	6✔
229	if len(value.shape) == 0:	6✔
230	mul, add = torch.tensor(1.0), torch.tensor(0.0)	6✔
231	if f"{key}_multiply" in scenario:	6✔
232	mul = scenario[f"{key}_multiply"]	6✔
233	if f"{key}_add" in scenario:	6✔
234	add = scenario[f"{key}_add"]	6✔
235
236	else:
237	add = torch.zeros_like(value)	6✔
238	mul = torch.ones_like(value)	6✔
239
240	if len(value.shape) == 1 and sec_vectors:	6✔
241	for s, ix in sec_vectors.items():	6✔
242	if f"{s}_{key}_multiply" in scenario:	6✔
243	mul = mul * (ix * scenario[f"{s}_{key}_multiply"])	6✔
244	if f"{s}_{key}_add" in scenario:	6✔
UNCOV 245	add = add + (ix * scenario[f"{s}_{key}_add"])	×
246
247	elif len(value.shape) == 2 and sec_matrices:	6✔
UNCOV 248	for rowcol, ix in sec_matrices.items():	×
UNCOV 249	s = f"{rowcol[0]}_{rowcol[1]}"	×
250
UNCOV 251	if f"{s}_{key}_multiply" in scenario:	×
252	mul = mul * (ix * scenario[f"{s}_{key}_multiply"])	×
UNCOV 253	if f"{s}_{key}_add" in scenario:	×
UNCOV 254	add = add + (ix * scenario[f"{s}_{key}_add"])	×
255
256	params[key] = value * mul + add	6✔
257
258	# Enforce adding-up constraints (differentiable). The resolver
259	# and constraint list are fetched fresh on every call so the
260	# step-time view of the parameter layout always matches the
261	# current Parameters instance.
262	constraints = self.parameters.get_constraints()	6✔
263	if constraints:	6✔
264	resolver = self.parameters.get_constraint_resolver()	6✔
265	for c in constraints:	6✔
266	c.apply(params, resolver)	6✔
267
268	return params	6✔
269
270	############################################################################
271	# Steady State
272	############################################################################
273
274	def compute_theoretical_steady_state(self, **kwargs):	6✔
275	"""Compute the theoretical steady state of the model.
276
277	This process generally follows the structure of the forward() function,
278	but instead of simulating the model, the steady state is computed at
279	each timestep. Therefore, (1) the model is initialized, and (2) for
280	each timestep the parameter and scenario information is passed to the
281	compute_theoretical_steady_state_per_step() function that computes the
282	steady state at that timestep.
283
284	Parameters
285	----------
286	**kwargs: dict
287	Additional keyword arguments.
288	"""
289	# Set the seed
290	torch.manual_seed(self.hyper["seed"])	6✔
291	self.numpy_rng = np.random.default_rng(self.hyper["seed"])	6✔
292
293	# Initialize the output tensors
294	self.state, self.history = self.variables.initialize_tensors()	6✔
295
296	# Initialize the model
297	info = f"(t=0...{self.hyper['timesteps_initialization']})"	6✔
298	logger.debug(f"Initializing model {info}")	6✔
299	self.initialize()	6✔
300
301	for t in range(self.hyper["timesteps_initialization"] + 1):	6✔
302	self.variables.record_state(t, self.state)	6✔
303
304	for t in range(self.hyper["timesteps_initialization"] + 1):	6✔
305	self.history = self.variables.update_history(self.state)	6✔
306
307	# Initialize prior from the initialization state
308	self.prior = self.state	6✔
309
310	# Compute the steady state
311	info = f"(t={self.hyper['timesteps_initialization'] + 1}...{self.hyper['timesteps']})"	6✔
312	logger.debug(f"Computing theoretical steady state {info}")	6✔
313
314	for t in range(	6✔
315	self.hyper["timesteps_initialization"] + 1, self.hyper["timesteps"]
316	):
317	self.state = self.variables.new_state()	6✔
318
319	# Get scenario series for this point in time
320	idx = torch.where(	6✔
321	torch.arange(self.hyper["timesteps"]) == t,
322	torch.ones(1),
323	torch.zeros(1),
324	)
325	scenario = {k: idx @ v for k, v in self.scenarios.items()}	6✔
326
327	# Apply parameter shocks
328	params = self.apply_parameter_shocks(t, scenario)	6✔
329
330	# Compute the steady state
331	self.compute_theoretical_steady_state_per_step(	6✔
332	t=t, params=params, scenario=scenario
333	)
334
335	# Store the outputs
336	self.variables.record_state(t, self.state)	6✔
337	self.history = self.variables.update_history(self.state)	6✔
338	self.prior = self.state	6✔
339
340	return None	6✔
341
342	def compute_theoretical_steady_state_per_step(self, **kwargs):	6✔
343	"""Compute the theoretical steady state of the model per step."""
344	raise NotImplementedError(
345	"Behavior.compute_theoretical_steady_state_per_step() to be implemented by model"
346	)
347
348	############################################################################
349	# Some Differentiable PyTorch Alternatives
350	############################################################################
351
352	def diffwhere(self, condition, x1, x2):	6✔
353	"""Where condition that is differentiable with respect to the condition.
354
355	Requires:
356	self.hyper['diffwhere'] = True
357	self.hyper['sigmoid_constant'] as a large number
358
359	Note: For non-NaN/inf, where(x > eps, z, y) is (x - eps > 0) * (z - y) + y,
360	so we can use the sigmoid function to approximate the where function.
361
362	Parameters
363	----------
364	condition : torch.Tensor
365	Condition to be evaluated expressed as x - eps
366	x1 : torch.Tensor
367	Value to be returned if condition is True
368	x2 : torch.Tensor
369	Value to be returned if condition is False
370	"""
371	sig = torch.sigmoid(torch.mul(condition, self.hyper["sigmoid_constant"]))	6✔
372	return torch.add(torch.mul(sig, torch.sub(x1, x2)), x2)	6✔
373
374	def tanhmask(self, x):	6✔
375	"""Convert a variable into 0 (x<0) and 1 (x>0)
376
377	Requires:
378	self.hyper['tanh_constant'] as a large number
379
380	Parameters
381	----------
382	x: torch.Tensor
383	The variable to be converted.
384
385	"""
386	kwg = {"dtype": torch.float64, "requires_grad": True}	6✔
387	return torch.div(	6✔
388	torch.add(
389	torch.ones(x.size(), **kwg),
390	torch.tanh(torch.mul(x, self.hyper["tanh_constant"])),
391	),
392	torch.tensor(2.0, **kwg),
393	)
394
395	def diffmin(self, x1, x2):	6✔
396	"""Smooth approximation to the minimum
397	B: https://mathoverflow.net/questions/35191/a-differentiable-approximation-to-the-minimum-function
398
399	Requires:
400	self.hyper['min_constant'] as a large number
401
402	Parameters
403	----------
404	x1: torch.Tensor
405	The first variable to be compared.
406	x2: torch.Tensor
407	The second variable to be compared.
408	"""
409	r = self.hyper["min_constant"]	6✔
410	pt1 = torch.exp(torch.mul(x1, -1 * r))	6✔
411	pt2 = torch.exp(torch.mul(x2, -1 * r))	6✔
412	return torch.mul(-1 / r, torch.log(torch.add(pt1, pt2)))	6✔
413
414	def diffmax(self, x1, x2):	6✔
415	"""Smooth approximation to the minimum
416	B: https://mathoverflow.net/questions/35191/a-differentiable-approximation-to-the-minimum-function
417
418	Requires:
419	self.hyper['max_constant'] as a large number
420
421	Parameters
422	----------
423	x1: torch.Tensor
424	The first variable to be compared.
425	x2: torch.Tensor
426	The second variable to be compared.
427	"""
428	r = self.hyper["max_constant"]	6✔
429	pt1 = torch.exp(torch.mul(x1, r))	6✔
430	pt2 = torch.exp(torch.mul(x2, r))	6✔
431	return torch.mul(1 / r, torch.log(torch.add(pt1, pt2)))	6✔
432
433	def diffmin_v(self, x):	6✔
434	"""Smooth approximation to the minimum. See diffmin
435
436	Parameters
437	----------
438	x: torch.Tensor
439	The variable to be converted.
440
441	Requires:
442	self.hyper['min_constant'] as a large number
443	"""
444	r = self.hyper["min_constant"]	6✔
445	temp = torch.exp(torch.mul(x, -1 * r))	6✔
446	return torch.mul(-1 / r, torch.log(torch.sum(temp)))	6✔
447
448	def diffmax_v(self, x):	6✔
449	"""Smooth approximation to the maximum for a tensor. See diffmax
450
451	Requires:
452	self.hyper['max_constant'] as a large number
453
454	Parameters
455	----------
456	x: torch.Tensor
457	The variable to be converted.
458	"""
459	r = self.hyper["max_constant"]	6✔
460	temp = torch.exp(torch.mul(x, r))	6✔
461	return torch.mul(1 / r, torch.log(torch.sum(temp)))	6✔
462
463
464	if __name__ == "__main__":
465	pass

KarlNaumann / MacroStat / 26623145945

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous