1.9. Optimal (optimized) simple rules

An optimal simple rule (OSR) keeps the model as written – including a simple, parametric feedback rule such as a Taylor rule – and chooses the rule coefficients that minimize a policy loss function. Unlike fully optimal policy (see Optimal policy), the model still has as many equations as endogenous variables: no equation is dropped, only some of its parameters are picked optimally.

1.9.1. The loss function

The loss is supplied to optimal_simple_rule at solve time, as a string in terms of the model’s variables and parameters (it is the planner’s objective, i.e. RISE maximises it – write a welfare measure, or minus a quadratic loss). A discount factor can be attached; it does not affect the ranking of unconditional losses, so it is often left at the default 0.99:

% minus a weighted quadratic loss in inflation, the output
% gap and the change in the policy rate
Loss = '-.5*(1*PAI^2 + .3*Y^2 + 0.9*DI^2)';

% with an explicit discount factor:
Loss = {0.99, '-.5*(1*PAI^2 + .3*Y^2 + 0.9*DI^2)'};

Important

The OSR objective is always supplied to optimal_simple_rule at solve time. It cannot be declared inside the model file: an @optimization_problem block in the .rs file is reserved strictly for optimal policy (the planner’s first-order conditions replace policy equations), which is a different solve concept from OSR (the policy equations are kept and only their coefficients are chosen optimally). See Optimal policy for the model-file block; here the loss is just a string argument.

The rule coefficients to optimize are listed the way estimated parameters are – through a priors structure giving each one a support (start value, lower bound, upper bound):

priors = struct();
priors.gam_lag = {0.60, 0.00,  1};     % {start, lb, ub}
priors.gam_y   = {0.50, 0.00, 10};

1.9.2. (Approximate) analytical loss

If the third argument is empty, optimal_simple_rule minimises the theoretical (unconditional) loss – computed from the model’s (perturbation-implied) second moments, with no simulation:

m    = dsge_model('Canonical_osr');
m    = set(m, parameters = p);

mest = optimal_simple_rule(m, Loss, [], estim_priors = priors);

print_estimation_results(mest)

1.9.3. Simulated loss

Passing a shocks database (or a struct / function handle producing one) as the third argument makes optimal_simple_rule minimize the conditional (simulated) loss along those shock paths instead – useful for nonlinear models or when you care about a specific scenario:

mest = optimal_simple_rule(m, Loss, shocks_db, estim_priors = priors);

You can evaluate the loss of any parameterization directly with calculate_loss(m, Loss) (unconditional) or calculate_loss(m, Loss, shocks_db) (conditional) – handy for comparing a hand-chosen rule with the optimized one.

calculate_loss returns the present value \(v(x,r) = x'V_r x + d_r\) of the planner’s objective. Because \(W\) is the Hessian of that objective, \(x'Wx = 2u\) for a quadratic loss, and for a stationary single-regime model the returned number satisfies exactly

\[|\mathrm{welf}| \;=\; \frac{2\,E[u]}{1-\beta}, \qquad E[u] = \tfrac12\sum_k \lambda_k \mathrm{Var}(x_k)\]

with the variances taken from theoretical_autocovariances. That identity is a convenient sanity check on any loss you write, and is asserted in the test suite (calculate_loss_identity).

The state \(x\) is measured in the coordinates of the model’s solution: deviations from the steady state for ordinary variables, and log deviations \(\log X_t - \log \bar X\) for variables declared @endogenous(log). The conditional loss converts the simulated levels accordingly, so it does not depend on how a log variable is normalized: rescaling a price level leaves it unchanged.

1.9.4. Switching models: how regimes are aggregated

In a Markov-switching model the unconditional loss is regime-specific: calculate_loss returns a 1 x h vector, one entry per regime.

optimal_simple_rule collapses that vector by summing the absolute values – every regime carries the same weight, regardless of how often it actually occurs. For a model whose bad regime is rare that is rarely what you want: a regime realised a fifth of the time counts as much as the one realised four fifths of it.

If you want a different aggregation – ergodic weights, or the loss conditional on a single regime – build the criterion yourself and pass it to indirect_inference, which takes a function handle (see OSR with a user-defined objective below). The ergodic distribution is available from the transition matrix:

Q  = m.model_data.state_space{1}.transition_matrices.Q;
pi = rise.engine.filtering.initial_markov_distribution(Q,true);

crit = pi(:).' * abs(welf(:));   % ergodic-weighted instead of summed

1.9.5. OSR and indirect inference

Because OSR is solved exactly like an estimation problem – an objective minimized over a set of parameters with bounds – it composes with RISE’s estimation machinery: the same optimizers (see Estimation) are available, and OSR can be combined with indirect inference (choosing rule coefficients to match empirical moments or impulse responses) by supplying the corresponding objective.

1.9.6. OSR with a user-defined objective

optimal_simple_rule takes the loss as an expression – the string (or {discount, string} cell) described above; it does not accept a function handle. If the objective you want cannot be written that way – for example matching empirical moments or impulse responses – use indirect inference directly, treating the rule coefficients as the estimated parameters. There the objective is a function handle, with the signature:

[critmin, retcode] = myobjective(m)

where m is the solved model object, critmin is the scalar criterion to be minimized, and retcode is 0 when the evaluation succeeded.

1.9.7. Policy frontiers

To trace a policy frontier – how the model’s volatilities trade off as a preference parameter (a weight in the loss, the persistence of the rule, …) is varied – use frontier:

% vary the weight on the output gap over a 50-point grid
% between 0 and 2
f = frontier(mest, 'lamb_y', [0, 2]);

% or pass an explicit grid
f = frontier(mest, 'lamb_y', linspace(0, 2, 21));

% use simulated rather than theoretical moments
f = frontier(mest, 'lamb_y', [0, 2], true);

f is a structure of the standard deviations of every model variable at each grid point (plus a stats__ sub-structure with the grid, the number of simulation periods, and a per-point return code); plot the columns of interest against each other to get the frontier.