This is the program file for Christopher Neely's article "Using Implied Volatility to
Measure Uncertainty About Interest Rates" in the May/June 2005 issue of
the Federal Reserve Bank of St. Louis Review.



****************************************************************************
****************************************************************************
****************************************************************************
****************************************************************************
/*
This program creates Table 1 for the article:

Christopher Neely, "Using Implied Volatility to Measure Uncertainty
About Interest Rates" Federal Reserve Bank of Saint Louis Review,
87(3), May/June 2005.

Programmer: Christopher Neely
Software: Gauss

This program creates summary statistics for on log futures price
changes (percent), annualized absolute log futures price changes, and
annualized IV and RV until expiry.  The rows show the total number of
observations in the sample, the non-missing observations, the mean,
the standard deviation, the maximum, the minimum, and the first 5
autocorrelations.  The standard error of the autocorrelations is
about 1/sqrt(T) ~ 0.016.

The futures data used are 3-month Eurodollar futures.  The data were
purchased from the Institute for Financial Markets (IFM).

*/

new;
library pgraph,arima,cml;
clear ufn,hu,grad2;
#include d:\CJN\Gaussprc\option.lib;
#include d:\CJN\Gaussprc\plotlab.src;

outwidth 256;

/* YYYYMMDD EXPIRY BDAYST_T USIRATE IV SETTLEFP ANNVOLTX LNDELTA */
inname1 = "d:\\cjn\\frb\\review\\implie~1\\data\\ED";
{d1} = OpenGaussSet(inname1);
{mat,decdate} = ymdtomatrix(d1[.,1]-19000000);

LnD = 100*d1[.,8];
LnD = missex(lnD, d1[.,3] .> lagn(d1[.,3],1));
LnD2 = 249*sqrt(LnD^2);
data = LnD~LnD2~100*d1[.,5|7];           @ DeltaLnFP DeltaLnFPSQ IV and RV @

TotalObs = ones(cols(data),1)*rows(data);
nobs = sumc(missrv(data+.0000001,0)./missrv(data+.0000001,1));
mean = sumc((missrv(data,0))./sumc(missrv(data,-999)./=-999)');
stddev = sqrt(sumc((missrv(data^2,0))./sumc(missrv(data^2,-999)./=-999)') - mean^2);

i= 1;
do while i <= cols(data);
if i == 1; acfX = acf(packr(data[.,i]),5,0); else; acfX = acfX~acf(packr(data[.,i]),5,0); endif;
i = i + 1;
endo;

stats = TotalObs~nobs~mean~stddev~maxc(data)~minc(data)~acfx';
coefnams = "TotalObs"|"Nobs"|"mean"|"stddev"|"max"|"min"|
            "AC1"|"AC2"|"AC3"|"AC4"|"AC5";

$coefnams~ftocv(stats',3,3);

outname = "D:\\cjn\\FRB\\REVIEW\\Implie~1\\out\\sumstats";
Let topnames = Stat LnDelta AnAbsLnD IVtilX RVtilX;
{z} = exportXLSandGauss(coefnams~stats',outname,topnames,1); 
z;






****************************************************************************
****************************************************************************
****************************************************************************
****************************************************************************
/*
This program creates Figures 2, 4, & 5, and Table 2 for the article:

Christopher Neely, "Using Implied Volatility to Measure Uncertainty
About Interest Rates" Federal Reserve Bank of Saint Louis Review,
87(3), May/June 2005.

Programmer: Christopher Neely
Software: Gauss

This program takes IV & RV for 3-month Eurodollar futures and:

_Figure 2_
Shows the empirical distributions of IV and changes in IV on 3-month
Eurodollar futures prices.

_Figure 4_
Displays 3-month Eurodollar IV and RV from 3/20/85 - 6/29/01.

_Figure 5_
Shows a scatter plot of (IV,RV) pairs along with the OLS fitted values
from Table 2.

_Table 2_
Shows the results of predicting 3-month Eurodollar RV with IV, as in
equation (8).

*/

new;
library pgraph,arima,cml;
clear ufn,hu,grad2;
#include d:\CJN\Gaussprc\option.lib;
#include d:\CJN\Gaussprc\plotlab.src;

outwidth 256;

/* YYYYMMDD EXPIRY BDAYST_T USIRATE IV SETTLEFP ANNVOLTX LNDELTA */
inname1 = "d:\\cjn\\frb\\review\\implie~1\\data\\ED";
{d1} = OpenGaussSet(inname1);
{mat,decdate} = ymdtomatrix(d1[.,1]-19000000);

datecv = miss(zeros(rows(d1),1),0); 

t = 1;
do while t <= rows(d1); 
datecv[t,1] = datestr( (1900+mat[t,1])|mat[t,2]|mat[t,3]|0); 
t = t + 1;
endo;

/*--- Figure 4 ---*/
graphset;
_pcolor = 8; _pmcolor = 8; _plwidth = 0|9; _pdate = 0;
_pltype = 1|6|2|3|4; _pcolor = 4|8|3|4|5; _paxht = .2; _pnumht = .2;
_pdate = 0; _pbox = 0; _pframe = 0; 
_ptitlht = .25;
_pframe = 0;
_plegctl = { 2 5 4 4 };
_plegstr = "Implied Volatility\000Realized Volatility";
@Title("Implied and Realized Volatility");@
ylabel("Percent");
xlabel("Date");
outgrf = "d:\\CJN\\FRB\\Review\\Implie~1\\grf\\Figure4.eps";
cmdstr = " -c=1" $+ " -cf=" $+ outgrf $+ " -w=5";
graphprt(cmdstr);
xy(decdate,100*d1[.,5|7]);  
/*----------------*/

pause(5); 

 graphset;
_pcolor = 8; _pmcolor = 8; _plwidth = 0|9; _pdate = 0;
_pltype = 1|6|2|3|4; _pcolor = 4|8|3|4|5; _paxht = .2; _pnumht = .2;
_ptitlht = .25;
_pframe = 0;
_plctrl = -1;
xlabel("Change in Futures Price");
ylabel("Change in IV");
outgrf = "d:\\CJN\\FRB\\Review\\Implie~1\\grf\\IV_DF.wmf";
ytics(-.4,.4,.2,0);
xtics(-.4,.6,.2,0);
xydata = ((d1[.,6]-lagn(d1[.,6],1))./100)~100*(d1[.,5]-lagn(d1[.,5],1));
packdat = packr( datecv~xydata~(d1[.,3] .> lagn(d1[.,3],1)) );
packdat = delif(packdat,packdat[.,4] .== 1);
datecv = packdat[.,1];
xydata = packdat[.,2:3];
plotlab(datecv, xydata, "Change in Futures Price", "Change in IV", "", 2, 2,outgrf);

xydata = delif(xydata,xydata[.,1] .== maxc(xydata[.,1]));
{vnam,m,b,stb,vc,stderr,sigma,cx,rsq,resid,dwstat}=ols(0,xydata[.,2],abs(xydata[.,1])); 

pause(5); 

/*--- Figure 2 ---*/
graphset;
begwind;
window(2,1,0);
setwind(1);
_pcolor = 8; _pmcolor = 8; _plwidth = 0|9; _pdate = 0;
_pltype = 1|6|2|3|4; _pcolor = 4|8|3|4|5; _paxht = .3; _pnumht = .3;
_ptitlht = .35;
_pframe = 0;
xlabel("Implied Volatility");
ylabel("Percent");
//xlabel("Date");
histp(100*d1[.,5],100); 
nextwind;
xlabel("Change in Implied Volatility");
histp(100*(d1[.,5]- lagn(d1[.,5],1)),100); 
outgrf = "d:\\CJN\\FRB\\Review\\Implie~1\\grf\\Figure2.eps";
cmdstr = " -c=1" $+ " -cf=" $+ outgrf $+ " -w=5";
graphprt(cmdstr);
endwind;   
/*----------------*/

/*  Load FOMC meeting dates */
FOMCname = "d:\\CJN\\FRB\\Review\\Implie~1\\data\\FOMC Meeting Dates Double.txt";
load fomc[]= ^FOMCname;
FOMC = reshape(fomc,rows(fomc)/3,3); 

/*  Load FF target changes */
FFname = "d:\\CJN\\FRB\\Review\\Implie~1\\data\\FFtargs.txt";
load FF[]= ^FFname;
FF = reshape(FF,rows(FF)/2,2); 
DFF = FF[2:rows(FF),2]-FF[1:rows(FF)-1,2];
DFF = delif(FF[2:rows(FF),1]~DFF,DFF .== 0);

/*  When are the largest changes in IV? */
Nchanges = 20;
IV = 100*d1[.,5];
RV = 100*d1[.,7];
DeltaIV = IV - lagn(IV,1);
DeltaRV = RV - lagn(RV,1);

SortDeltaIV = rev(sortc(abs(DeltaIV)~seqa(1,1,rows(DeltaIV)),1));
SortDeltaRV = rev(sortc(abs(DeltaRV)~seqa(1,1,rows(DeltaRV)),1));

SortDeltaIV = DeltaIV[SortDeltaIV[.,2],1]~d1[SortDeltaIV[.,2],1];

SortIVdates = SortDeltaIV[.,2]-10000*trunc(SortDeltaIV[.,2]./10000);
bad = maxc((sortIVdates .== 301~601~901~1201)') .== 1;
SortDeltaIV = delif(SortDeltaIV,bad);

fomcindx = indnv(SortDeltaIV[1:Nchanges,2],DFF[.,1]|DFF[.,1]+1);
deltaFindx = indnv(SortDeltaIV[1:Nchanges,2],D1[.,1]);

SortDeltaRV = DeltaRV[SortDeltaRV[.,2],1]~d1[SortDeltaRV[.,2],1];

SortRVdates = SortDeltaRV[.,2]-10000*trunc(SortDeltaRV[.,2]./10000);
bad = maxc((sortRVdates .== 301~601~901~1201)') .== 1;
SortDeltaRV = delif(SortDeltaRV,bad);

d1[deltaFindx[1:20,.],8];       @ Futures price changes for largest changes in IV @
pv = meanc(packr(abs(d1[.,8])) .> abs(d1[deltaFindx[1:20,.],8])');  @ P-value @

/*----------------- Big News ----------------------*/
outname = "d:\\CJN\\FRB\\Review\\Implie~1\\out\\bignews.txt";
output file = ^outname reset;
"Largest changes in IV (size & dates), FOMCINDX & RV (size & dates), deltaF & PV";
SortDeltaIV[1:Nchanges,.]~FOMCINDX~SortDeltaRV[1:Nchanges,.]~d1[deltaFindx[1:20,.],8]~pv;
output off; 


/*----------------- Plot just the largest changes ----------------*/
smallindx = indnv(SortDeltaIV[1:Nchanges,2],d1[.,1]);
smallindx = smallindx[1:10,.];
xydata = ((d1[.,6]-lagn(d1[.,6],1))./100)~100*(d1[.,5]-lagn(d1[.,5],1));
packdat = packr(datecv[smallindx,.]~xydata[smallindx,.]);
datecv = packdat[.,1];
xydata = packdat[.,2:3];


graphset;
_pcolor = 8; _pmcolor = 8; _plwidth = 0|9; _pdate = 0;
_pltype = 1|6|2|3|4; _pcolor = 4|8|3|4|5; _paxht = .2; _pnumht = .2;
_ptitlht = .25;
_pframe = 0;
_plctrl = -1;
xlabel("Change in Futures Price");
ylabel("Change in IV");
outgrf = "d:\\CJN\\FRB\\Review\\Implie~1\\grf\\IV_DF_small.wmf";

/*  cmdstr = " -c=8" $+ " -cf=" $+ outgrf $+ " -w=15";
graphprt(cmdstr);
xy(100*(d1[.,6]-lagn(d1[.,6],1)),100*(d1[.,5]-lagn(d1[.,5],1)));   */

ytics(-.6,.6,.6,0);
xtics(-.5,1.2,.8,0);

plotlab(datecv, xydata, "Change in Futures Price", "Change in IV", "", 2, 2,outgrf);

/*-------------- Predict RV with IV (Table 2) --------------*/
outname = "d:\\CJN\\FRB\\Review\\Implie~1\\out\\OLSRVonIV.txt";
_olsres = 1;

olsdata = packr(d1[.,3|5|7]);       @ BusT_t, IV, RV @
ydat = 100*olsdata[.,3];            @ Realized Volatility @
xdat = 100*olsdata[.,2];            @ IV @
timetoexpiry = olsdata[.,1]; 

output file = ^outname reset;
{vnam,m,b1,stb,vc,stderr,sigma,cx,rsq,resid,dw}=ols(0,ydat,xdat); 

{HCVC1,VC1} = JorionSE(packr(ones(rows(xdat),1)~xdat),resid,timetoexpiry);
stderr1 = sqrt(diag(HCVC1));
BigR = eye(2);  littler = 0|1;
Wald1 = (BigR*b1-littler)'inv(BigR'HCVC1*BigR)*(BigR*b1-littler);
WaldPv = cdfchic(Wald1,2);
output off; 

LHS = "a"|"(s.e.)"|"B"|"(s.e.)"|"Wald"|"WaldPV"|"Nobs"|"RSQ";
outdata = b1[1]|stderr1[1]|b1[2]|stderr1[2]|wald1|waldpv|rows(ydat)|RSQ;
outdata = LHS~outdata;

outname = "d:\\CJN\\FRB\\Review\\Implie~1\\out\\RVonIV";
let topnames = var stat;
{z} = exportXLSandGauss(outdata,outname,topnames,1);


/*------- Plot the X,Y coordinates RV vs IV and fitted line (Figure 5) -------*/
fitline = b1[1] + b1[2]*xdat;
line45 = seqa(0,maxc(ydat)/rows(xdat),rows(xdat));

graphset;
_pcolor = 8; _pmcolor = 8; _plwidth = 0|9|9; _pdate = 0;
_pltype = 1|6|2|3|4; _pcolor = 4|8|3|4|5; _paxht = .2; _pnumht = .2;
_ptitlht = .25;
_pframe = 0;
_plctrl = {-1 0 0};
_psym = meanc(xdat)~meanc(ydat)~4~12.0~8~1~9;

/* Mx7 matrix, M extra symbols will be plotted. 
 [M,1] x location.  [M,2] y location.  [M,3] symbol type. (See _pstype, earlier.) 
 [M,4] symbol height. if this is 0, a default height of 5.0 will be used. 
 [M,5] symbol color, see "Colors," page 22-12.  [M,6] type of coordinates: 1 plot 2 inch  
 [M,7] line thickness. A value of zero is normal line width.   */

ylabel("Realized Volatility until Expiry");
xlabel("Implied Volatility");
outgrf = "d:\\CJN\\FRB\\Review\\Implie~1\\grf\\Figure5.eps";
cmdstr = " -c=1" $+ " -cf=" $+ outgrf $+ " -w=5";
graphprt(cmdstr);
xy(xdat~xdat~line45,ydat~fitline~line45);  


****************************************************************************
****************************************************************************
****************************************************************************
****************************************************************************
/*
This program creates figure 3 for the article:

Christopher Neely, "Using Implied Volatility to Measure Uncertainty
About Interest Rates" Federal Reserve Bank of Saint Louis Review,
87(3), May/June 2005.

Programmer: Christopher Neely
Software: Gauss

This program takes 3 month Eurodollar rates (from the Federal Reserve
Board's H.15 release), and the Fed Funds target rate, and plots them
from 1/1/84 - 7/25/03.

*/

new;
library pgraph,arima,cml;
clear ufn,hu,grad2;
#include d:\CJN\Gaussprc\option.lib;

inname1 = "d:\\cjn\\frb\\review\\implie~1\\data\\ED3mobidH15.txt";
load m1[] = ^inname1;
m1 = reshape(m1,rows(m1)/2,2);
m1 = delif(m1,m1[.,1] .< 19840000);
YYYYMMDD1 = m1[.,1];
ED3mo = missex(m1[.,2],m1[.,2] .== 0);
{mat1,decdate1} = ymdtomatrix(yyyymmdd1);


inname2 = "d:\\cjn\\frb\\review\\implie~1\\data\\USFFtarg.txt";
load m2[] = ^inname2;
m2 = reshape(m2,rows(m2)/2,2);
YYYYMMDD2 = m2[.,1];
fftargets = m2[.,2];
{mat2,decdate2} = ymdtomatrix(yyyymmdd2);

xdata = miss(zeros(maxc(rows(decdate1)|rows(decdate2)),2),0);
ydata = xdata;

xdata[1:rows(decdate1),1] = decdate1;
xdata[1:rows(decdate2),2] = decdate2;
ydata[1:rows(decdate1),1] = ED3mo;
ydata[1:rows(decdate2),2] = FFtargets;

graphset;
_pcolor = 8; _pmcolor = 8; _plwidth = 0|9; _pltype = 1|6; _pcolor = 4|8;
_pdate = 0; _pbox = 0; _pframe = 0; 
xtics(1984,2004,5,0);
ylabel("Federal Funds Targets and 3-month Eurodollar");
xlabel("Date");
outgrf = "d:\\CJN\\FRB\\Review\\Implie~1\\grf\\Figure3.eps";
cmdstr = " -c=1" $+ " -cf=" $+ outgrf $+ " -w=5";
graphprt(cmdstr);
xy(xdata,ydata);