diff --git a/HARK/ConsumptionSaving/ConsPortfolioFrameModel.py b/HARK/ConsumptionSaving/ConsPortfolioFrameModel.py new file mode 100644 index 000000000..44b483598 --- /dev/null +++ b/HARK/ConsumptionSaving/ConsPortfolioFrameModel.py @@ -0,0 +1,147 @@ +""" +This file contains classes and functions for representing, +solving, and simulating agents who must allocate their resources +among consumption, saving in a risk-free asset (with a low return), +and saving in a risky asset (with higher average return). + +This file also demonstrates a "frame" model architecture. +""" +import numpy as np +from scipy.optimize import minimize_scalar +from copy import deepcopy +from HARK import HARKobject, NullFunc, Frame, FrameAgentType # Basic HARK features +from HARK.ConsumptionSaving.ConsIndShockModel import ( + IndShockConsumerType, # PortfolioConsumerType inherits from it + ValueFunc, # For representing 1D value function + MargValueFunc, # For representing 1D marginal value function + utility, # CRRA utility function + utility_inv, # Inverse CRRA utility function + utilityP, # CRRA marginal utility function + utility_invP, # Derivative of inverse CRRA utility function + utilityP_inv, # Inverse CRRA marginal utility function + init_idiosyncratic_shocks, # Baseline dictionary to build on +) +from HARK.ConsumptionSaving.ConsGenIncProcessModel import ( + ValueFunc2D, # For representing 2D value function + MargValueFunc2D, # For representing 2D marginal value function +) + +from HARK.ConsumptionSaving.ConsPortfolioModel import ( + init_portfolio, + solveConsPortfolio, + PortfolioConsumerType, + PortfolioSolution +) + +from HARK.distribution import combineIndepDstns +from HARK.distribution import Lognormal, MeanOneLogNormal, Bernoulli # Random draws for simulating agents +from HARK.interpolation import ( + LinearInterp, # Piecewise linear interpolation + CubicInterp, # Piecewise cubic interpolation + LinearInterpOnInterp1D, # Interpolator over 1D interpolations + BilinearInterp, # 2D interpolator + ConstantFunction, # Interpolator-like class that returns constant value + IdentityFunction, # Interpolator-like class that returns one of its arguments +) + +class PortfolioConsumerFrameType(FrameAgentType, PortfolioConsumerType): + """ + A consumer type with a portfolio choice, using Frame architecture. + + A subclass of PortfolioConsumerType for now. + This is mainly to keep the _solver_ logic intact. + """ + + # values for aggregate variables + # to be set when simulation initializes. + # currently not doing anything because still using old + # initializeSim() + aggregate_init_values = { + 'PlvlAggNow' : 1.0 + } + + def birth_aNrmNow(self, N): + return Lognormal( + mu=self.aNrmInitMean, + sigma=self.aNrmInitStd, + seed=self.RNG.randint(0, 2 ** 31 - 1), + ).draw(N) + + def birth_pLvlNow(self, N): + pLvlInitMeanNow = self.pLvlInitMean + np.log( + self.state_now["PlvlAggNow"] + ) # Account for newer cohorts having higher permanent income + + return Lognormal( + pLvlInitMeanNow, + self.pLvlInitStd, + seed=self.RNG.randint(0, 2 ** 31 - 1) + ).draw(N) + + + # values to assign to agents at birth + birth_values = { + 'ShareNow' : 0, + 'AdjustNow' : False, + 'aNrmNow' : birth_aNrmNow, + 'pLvlNow' : birth_pLvlNow + } + + def transition_ShareNow(self, **context): + ShareNow = np.zeros(self.AgentCount) + np.nan + + # Loop over each period of the cycle, getting controls separately depending on "age" + for t in range(self.T_cycle): + these = t == self.t_cycle + + # Get controls for agents who *can* adjust their portfolio share + those = np.logical_and(these, self.shocks['AdjustNow']) + + ShareNow[those] = self.solution[t].ShareFuncAdj(self.state_now['mNrmNow'][those]) + + # Get Controls for agents who *can't* adjust their portfolio share + those = np.logical_and( + these, + np.logical_not(self.shocks['AdjustNow'])) + ShareNow[those] = self.solution[t].ShareFuncFxd( + context['mNrmNow'][those], ShareNow[those] + ) + + # redundant for now + self.controls["ShareNow"] = ShareNow + + return ShareNow + + def transition_cNrmNow(self, **context): + cNrmNow = np.zeros(self.AgentCount) + np.nan + ShareNow = self.controls["ShareNow"] + + # Loop over each period of the cycle, getting controls separately depending on "age" + for t in range(self.T_cycle): + these = t == self.t_cycle + + # Get controls for agents who *can* adjust their portfolio share + those = np.logical_and(these, context['AdjustNow']) + cNrmNow[those] = self.solution[t].cFuncAdj(context['mNrmNow'][those]) + + # Get Controls for agents who *can't* adjust their portfolio share + those = np.logical_and( + these, + np.logical_not(context['AdjustNow'])) + cNrmNow[those] = self.solution[t].cFuncFxd( + context['mNrmNow'][those], ShareNow[those] + ) + + # Store controls as attributes of self + # redundant for now + self.controls["cNrmNow"] = cNrmNow + + return cNrmNow + + frames = { + ('RiskyNow','AdjustNow') : PortfolioConsumerType.getShocks, + ('pLvlNow', 'PlvlAggNow', 'bNrmNow', 'mNrmNow') : PortfolioConsumerType.getStates, + ('ShareNow') : transition_ShareNow, + ('cNrmNow') : transition_cNrmNow, + ('aNrmNow', 'aNrmNow') : PortfolioConsumerType.getPostStates + } diff --git a/HARK/ConsumptionSaving/ConsPortfolioModel.py b/HARK/ConsumptionSaving/ConsPortfolioModel.py index f938777d0..b976886f7 100644 --- a/HARK/ConsumptionSaving/ConsPortfolioModel.py +++ b/HARK/ConsumptionSaving/ConsPortfolioModel.py @@ -444,10 +444,7 @@ def simBirth(self, which_agents): None """ IndShockConsumerType.simBirth(self, which_agents) - # Checking for control variable attribute here - # because we have not namespaced controls yet - if hasattr(self, 'ShareNow'): - self.ShareNow[which_agents] = 0 + self.controls['ShareNow'][which_agents] = 0 # here a shock is being used as a 'post state' self.shocks['AdjustNow'][which_agents] = False diff --git a/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py b/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py new file mode 100644 index 000000000..d5b699ce3 --- /dev/null +++ b/HARK/ConsumptionSaving/tests/test_ConsPortfolioFrameModel.py @@ -0,0 +1,120 @@ +import HARK.ConsumptionSaving.ConsPortfolioFrameModel as cpfm +import numpy as np +import unittest + + +class PortfolioConsumerTypeTestCase(unittest.TestCase): + def setUp(self): + # Create portfolio choice consumer type + self.pcct = cpfm.PortfolioConsumerFrameType() + self.pcct.cycles = 0 + + # Solve the model under the given parameters + + self.pcct.solve() + +class UnitsPortfolioConsumerTypeTestCase(PortfolioConsumerTypeTestCase): + def test_simOnePeriod(self): + + self.pcct.T_sim = 30 + self.pcct.AgentCount = 10 + self.pcct.track_vars += ['aNrmNow'] + self.pcct.initializeSim() + + self.assertFalse( + np.any(self.pcct.shocks['AdjustNow']) + ) + + self.pcct.simOnePeriod() + + self.assertAlmostEqual( + self.pcct.controls["ShareNow"][0], + 0.8627164488246847 + ) + self.assertAlmostEqual( + self.pcct.controls["cNrmNow"][0], + 1.67874799 + ) + +class SimulatePortfolioConsumerTypeTestCase(PortfolioConsumerTypeTestCase): + + def test_simulation(self): + + self.pcct.T_sim = 30 + self.pcct.AgentCount = 10 + self.pcct.track_vars += [ + 'mNrmNow', + 'cNrmNow', + 'ShareNow', + 'aNrmNow', + 'RiskyNow', + 'RportNow', + 'AdjustNow', + 'PermShkNow', + 'bNrmNow' + ] + self.pcct.initializeSim() + + self.pcct.simulate() + + self.assertAlmostEqual( + self.pcct.history['mNrmNow'][0][0], 9.70233892 + ) + + self.assertAlmostEqual( + self.pcct.history['cNrmNow'][0][0], 1.6787479894848298 + ) + + self.assertAlmostEqual( + self.pcct.history['ShareNow'][0][0], 0.8627164488246847 + ) + + self.assertAlmostEqual( + self.pcct.history['aNrmNow'][0][0], 8.023590930905383 + ) + + self.assertAlmostEqual( + self.pcct.history['AdjustNow'][0][0], 1.0 + ) + + + # the next period + self.assertAlmostEqual( + self.pcct.history['RiskyNow'][1][0], 0.8950304697526602 + ) + + self.assertAlmostEqual( + self.pcct.history['RportNow'][1][0], 0.9135595661654792 + ) + + self.assertAlmostEqual( + self.pcct.history['AdjustNow'][1][0], 1.0 + ) + + self.assertAlmostEqual( + self.pcct.history['PermShkNow'][1][0], 1.0050166461586711 + ) + + self.assertAlmostEqual( + self.pcct.history['bNrmNow'][1][0], 7.293439643953855 + ) + + self.assertAlmostEqual( + self.pcct.history['mNrmNow'][1][0], 8.287859049575047 + ) + + self.assertAlmostEqual( + self.pcct.history['cNrmNow'][1][0], 1.5773607434989751 + ) + + self.assertAlmostEqual( + self.pcct.history['ShareNow'][1][0], 0.9337608822146805 + ) + + self.assertAlmostEqual( + self.pcct.history['aNrmNow'][1][0], 6.710498306076072 + ) + + self.assertAlmostEqual( + self.pcct.history['aNrmNow'][15][0], 5.304746367434934 + ) diff --git a/HARK/core.py b/HARK/core.py index cb9e09dde..c1b7d8fb6 100644 --- a/HARK/core.py +++ b/HARK/core.py @@ -454,8 +454,14 @@ def initializeSim(self): if self.state_now[var] is None: self.state_now[var] = copy(blank_array) - #elif self.state_prev[var] is None: - # self.state_prev[var] = copy(blank_array) + for var in self.controls: + if self.controls[var] is None: + self.controls[var] = copy(blank_array) + + for var in self.shocks: + if self.shocks[var] is None: + self.shocks[var] = copy(blank_array) + self.t_age = np.zeros( self.AgentCount, dtype=int ) # Number of periods since agent entry @@ -806,6 +812,183 @@ def clearHistory(self): self.history[var_name] = np.empty((self.T_sim, self.AgentCount)) + np.nan +class Frame(): + """ + """ + + def __init__( + self, + target, + scope, + default = None, + transition = None, + objective = None + ): + """ + """ + + self.target = target # tuple of variables + self.scope = scope # tuple of variables + self.default = default # default value used in simBirth; a dict + self.transition = transition # for use in simulation + self.objective = objective # for use in solver + + +class FrameAgentType(AgentType): + """ + A variation of AgentType that uses Frames to organize + its simulation steps. + + Frames allow for state, control, and shock resolutions + in a specified order, rather than assuming that they + are resolved as shocks -> states -> controls -> poststates. + + Attributes + ---------- + state_vars : list of string + The string labels for this AgentType's model state variables. + """ + + # frames property + frames = [ + Frame( + ('y'),('x'), + transition = lambda x: x^2 + ) + ] + + def simOnePeriod(self): + """ + Simulates one period for this type. Calls the methods getMortality(), getShocks() or + readShocks, getStates(), getControls(), and getPostStates(). These should be defined for + AgentType subclasses, except getMortality (define its components simDeath and simBirth + instead) and readShocks. + + Parameters + ---------- + None + + Returns + ------- + None + """ + if not hasattr(self, "solution"): + raise Exception( + "Model instance does not have a solution stored. To simulate, it is necessary" + " to run the `solve()` method of the class first." + ) + + # Mortality adjusts the agent population + self.getMortality() # Replace some agents with "newborns" + + # state_{t-1} + for var in self.state_now: + self.state_prev[var] = self.state_now[var] + # note: this is not type checked for aggregate variables. + self.state_now[var] = np.empty(self.AgentCount) + + # transition the variables in the frame + for frame in self.frames: + self.transition(frame) + + # Advance time for all agents + self.t_age = self.t_age + 1 # Age all consumers by one period + self.t_cycle = self.t_cycle + 1 # Age all consumers within their cycle + self.t_cycle[ + self.t_cycle == self.T_cycle + ] = 0 # Resetting to zero for those who have reached the end + + def simBirth(self, which_agents): + """ + Makes new agents for the simulation. + Takes a boolean array as an input, indicating which + agent indices are to be "born". + + Populates model variable values with value from `init` + property + + Parameters + ---------- + which_agents : np.array(Bool) + Boolean array of size self.AgentCount indicating which agents should be "born". + + Returns + ------- + None + """ + N = np.sum(which_agents) + + for frame in self.frames: + for var in frame.target: + + if callable(frame.default[var]): + value = frame.default[var](self, N) + else: + value = frame.default[var] + + if var in self.state_now: + self.state_now[var][which_agents] = value + elif var in self.controls: + self.controls[var][which_agents] = value + elif var in self.shocks: + self.shocks[var][which_agents] = value + + # from ConsIndShockModel. Needed??? + self.t_age[which_agents] = 0 # How many periods since each agent was born + self.t_cycle[ + which_agents + ] = 0 # Which period of the cycle each agent is currently in + + def transition_frame(self, frame): + """ + Updates the model variables in `target` + using the `transition` function. + + The transition function will use current model + variable state as arguments. + """ + + # build a context object based on model state variables + # and 'self' reference for 'global' variables + context = {'self' : self} + context.update{self.shocks} + context.update{self.controls} + context.update{self.state_now} + + # a method for indicating that a 'previous' version + # of a variable is intended. + # Perhaps store this in a separate notation.py module + def decrement(var_name): + return var_name + '_' + + # use special notation for the 'previous state' variables + context.update({ + decrement(var) : state_prev[var] + for var + in state_prev + + }) + + # limit context to scope of frame + local_context = { + var : context[var] + for var + in frame.scope + } + + if frame.transition is not None: + new_values = frame.transition( + **local_context + ) + else: + raise Exception(f"Frame has None for transition: {frame}") + + # because the context was a shallow update, + # the model values can be modified directly(?) + for i in enumerate(frame_target): + context[target[i]] = new_values[i] + + def solveAgent(agent, verbose): """ Solve the dynamic model for one agent type