@
This file contains the Gauss programs that calculate inflation-output 
variablity for the case of New-Keynesian price-level and inflation targetting 
(figures 1-5)in "What do New-Keynesian Phillips Curves Imply for Price-Level 
Targeting?" by Robert Dittmar and William T. Gavin in the March/April 2000 
issue of the Federal Reserve Bank of St. Louis Review. This file contains two
 programs.  The first program calculates inflation-output variability in the 
case of price-level targeting. The second program calculates inflation-output 
variability in the case of inflation targeting. 
@

/*-----------------------------------------------------------------------------*/
/* Program to calculate inflation-output variability for price-level targeting */
/*-----------------------------------------------------------------------------*/
@ Parameters @
r = (1.04)^.25 - 1;       @ Interest rate                           @
xi = 1.0^.25;             @ Growth rate of labor prod.              @
beta = xi/(1+r);          @ Time preference                         @
mustar = xi*(1.0497^.25); @ Steady state growth rate of mon. supply @
pistar = mustar/xi;       @ Steady state inflation rate             @

@ Parameters from the monetary authorities optimization problem @
lambda = .5; @ Weight on output gap       @
rho = .9;    @ Phillip's curve parameters @
alpha = .5;

/* Calculation of coeffients on policy rule for y[t] */
PROC (1) = BCalc(A);

LOCAL B;

B = ZEROS(3,1);
B[1] = -alpha*(A[1] - 1)^2/(1 + alpha*A[2]);
B[2] = (rho - alpha*A[2]*(A[1] - 2))/(1 + alpha*A[2]);
B[3] = (1 - alpha*A[3]*(A[1] - 2))/(1 + alpha*A[2]);

RETP(B);

ENDP;


/* Coefficient matching conditions */
PROC (1) = coeffMtch(A);

LOCAL B,determ,lams,nuMat,matchEqs,i,alfs,condition,matrini;

B = BCalc(A);

matchEqs = {};
i = 1;
DO WHILE i <= 3;
  condition = B[i]*(A[2]^2*B[1]*beta^3*rho^2 + alpha*(-2 + A[1]*beta)*(-2 +
    B[2]*beta)*lambda*(-1 + A[1]*beta*rho) - A[2]*beta*(-1 + A[1]*beta*rho)*
    (-1 + alpha*B[1]*beta*lambda + B[2]*beta*rho)) + A[i]*(-2 -
    4*alpha*B[1]*beta*lambda*rho + A[2]*alpha*B[1]^2*beta^3*lambda*rho +
    2*A[1]*beta*(1 + alpha*B[1]*beta*lambda)*rho + 2*A[2]*B[1]*beta^2*rho^2 +
    B[2]^2*beta^2*rho*(-1 + A[1]*beta*rho) - B[2]*beta*(-1 +
    (-2 + A[1]*beta)*(1 + alpha*B[1]*beta*lambda)*rho + beta*(2*A[1] +
    A[2]*B[1]*beta)*rho^2));
  matchEqs = matchEqs|condition;
  i = i + 1;
ENDO;

RETP(matchEqs);

ENDP;

varList = {};
w = 1/100;
DO WHILE w <= 99/100;

lambda = w/(1 - w);

xk = ZEROS(3,1);
Ak = GRADP(&coeffMtch,xk);
Fk = coeffMtch(xk);

@ Globally convergent (?) Broyden's method @
_stepTol = 3*2.5*10^(-24); _stepTol;
solved = 0;
i = 1;
DO UNTIL solved;
  sk = -Fk/Ak;
  norm = MAXC(ABS(sk));
  @"Max. Broyden step" norm;@
  IF sk'sk <= _stepTol;

  solved = 1;

  ELSE;

  done = 0;
  DO UNTIL done;

    xkP1 = xk + sk;
    FkP1 = coeffMtch(xkP1);
    cond1 = FkP1'FkP1/2;
    cond2 = Fk'Fk/2 + .0001*Fk'*Ak*sk;

    IF ISMISS(FkP1'FkP1);
      "Infeasible coefficient values";
      lam = .1;
      sk = lam*sk;
    ELSEIF (cond1 > cond2);
      lam = -Fk'*Ak*sk/(FkP1'FkP1 - Fk'Fk - 2*Fk'*Ak*sk);
      lam = MAXC(lam|.1);
      lam = MINC(lam|.5);
      @"Backtracking" lam;
      backNum = backNum + 1;@
      sk = lam*sk;
    ELSE;
      done = 1;
    ENDIF;
  ENDO;

  xk = xkP1;
  Fk = FkP1;
  Ak = GRADP(&coeffMtch,xk);

  @ Relative gradient criteria @
  relGrad = 2*(Ak'*Fk).*xk/(Fk'Fk);

  @"Step #" i;
  "F(x)'F(x)/2" (Fk'Fk)/2;
  "Rel. grad." MAXC(ABS(relGrad));

  SAVE xk; SAVE Ak;@

  i = i + 1;

  ENDIF;

ENDO;

/* Coefficients on policy rules */
A = xk;
B = BCalc(A);
@eig(a[1:2]'|b[1:2]'); wait;@

unConVarMat = (1~-(A[1]-1)^2~-2*A[2]*(A[1]-1)~-A[2]^2)|
  (0~1-A[1]^2~-2*A[1]*A[2]~-A[2]^2)|
  (0~-A[1]*B[1]~1-(A[2]*B[1]+A[1]*B[2])~-A[2]*B[2])|
  (0~-B[1]^2~-2*B[1]*B[2]~1-B[2]^2);
RHS = (A[3]^2|A[3]^2|A[3]*B[3]|B[3]^2);
vars = RHS/unConVarMat;

@ Comparision to neo-Classical Phillips curve trade-offs @
Kap = SQRT((1 + rho)*vars[4]/2);
pLam = (1 - beta*rho*rho)*Kap/(alpha*(1 - alpha*Kap)); pLam;
IF pLam < 0;
  varYNC = MISS(0,0);
ELSE;
  varYNC = ((1 - beta*rho*rho)^2)/((1 - rho*rho)*(1 - beta*rho*rho +
    alpha*alpha*pLam)^2);
ENDIF;

varList = varList|(vars[4]~vars[1]);
@Column 1 of varList is inflation variability and column 2 is output variability@

w = w + 1/100;

ENDO;


/*
---------------------------------------------------------------------------------
---------------------------------------------------------------------------------
---------------------------------------------------------------------------------
*/


/*-----------------------------------------------------------------------------*/
/* Program to calculate inflation-output variability for inflation targeting   */
/*-----------------------------------------------------------------------------*/
@ Parameters @
r = (1.04)^.25 - 1;       @ Interest rate                           @
xi = 1.0^.25;             @ Growth rate of labor prod.              @
beta = xi/(1+r);          @ Time preference                         @
mustar = xi*(1.0497^.25); @ Steady state growth rate of mon. supply @
pistar = mustar/xi;       @ Steady state inflation rate             @

varList = {};
w = 1/100;
DO WHILE w <= 99/100;

lambda = w/(1 - w);

@ Parameters from the monetary authorities optimization problem @
@lambda = .5;@ @ Weight on output gap       @
rho = .9;    @ Phillip's curve parameters @
alpha = .5;

@ Calculation of decision rules coefficients @
a = alpha*(1 - beta*rho);
b = (1 - beta*rho*rho + alpha*alpha*lambda);
c = alpha*lambda*rho;
A2 = POLYROOT(a|b|c); a|b|c;
A2 = A2~((rho + alpha*A2)./(1 + alpha*A2));
A2 = SELIF(A2,ABS(A2[.,2]) .<= 1);

B2 = A2[.,2]; A2 = A2[.,1];
A3 = (1 + alpha*A2)/(beta*rho*A2 - alpha*lambda) - alpha;
A3 = 1/A3;
B3 = (1 + alpha*A3)/(1 + alpha*A2);

@ Unconditional variances of output and inflation @
varY = B3*B3/(1 - B2*B2);
varPi = A2*A2*varY + A3*A3;

@ Comparision to neo-Classical Phillips curve trade-offs @
Kap = SQRT((1 - rho*rho)*varPi);
pLam = (1 - beta*rho*rho)*Kap/(alpha*(1 - alpha*Kap));
varYNC = ((1 - beta*rho*rho)^2)/((1 - rho*rho)*(1 - beta*rho*rho +
  alpha*alpha*pLam)^2);

varList = varList|(varPi~varY);
@Column 1 of varList is inflation variability and column 2 is output variability@

((A2~B2)|(A3~B3));

w = w + 1/100;

ENDO;