Import & Analyze Natural Gas Futures Options & Historical Prices

This is the first script of 4 in the Natural Gas Storage Valuation case study. In this script, we import a data set of historical futures prices for natural gas, pre-process them and compute return statistics. We also compute implied volatilities for a set of current NG futures option prices. Instructions for obtaining the data set are provided in the Readme document.


Define Parameters

Set the valuation date and start date for historical statistics, and discounting rate for option pricing

valDate = '2014-05-15';
histStartDate = '2004-01-01';
discRate = .03;
addpath Utilities

Import Historical Forward Prices

Import "tall" table of settle prices of historical natural gas futures contracts. Each row in this table corresponds to one historical observation of a single futures price.

data = readtable('Data/NGFuturesPriceHistory.xlsx');

Pre-process Futures Prices

Filter the data set to keep prices between the start date and valuation date. Then, use unstack to reshape the data set so that each row contains prices for the full futures curve on that date.

data.Date = datenum(data.Date, getLocalDateFormat);
data = data(data.Date >= datenum(histStartDate) & ...
            data.Date <= datenum(valDate), :);

% "Pivot" settle prices on futures symbol (tall -> wide)
curveSet = unstack(data, 'Settle', 'Symbol');
expDate = futSymbolToExp(curveSet.Properties.VariableNames(2:end));
[obsDate, idx1] = sort(curveSet.Date);
[expDate, idx2] = sort(expDate);
curveSet = curveSet(idx1,[1;1+idx2]);

% Extract prices into a matrix
fwdPrice = table2array(curveSet(:,2:end));

Visualize Historical Futures Curves

Use a surface plot to visualize historical curves and detect missing data

plotForwardSeries(obsDate, expDate, fwdPrice);

Remove long dated contracts

For the purpose of our analysis, we may not need more than one or two years worth of contracts. The others may be removed.

fwdPrice = keepExpirations(fwdPrice, 1:24,[],'removeEmptyCols');
expDate = expDate(1:size(fwdPrice,2));
%plotForwardSeries(obsDate, expDate, fwdPrice);

save Data/FwdCurveHistory fwdPrice expDate obsDate

Compute Promptness (Maturity) Returns

Compute log returns for each contract. Then, convert from "Contract" (columns representing contracts) to "Maturity" (columns representing promptness) representation. This involves shifting the curve every time a contract expires.

f = log(fwdPrice);
df = diff(f, 1, 1); % First difference along first dimension
df = collapseFwdMatrix(df); % Contract to promptness

plotForwardSeries(obsDate(2:end), 1:size(df,2), df);

Compute Statistics of Forward Curve Return Series

Compute means, variance, correlation and covariance for the forward curve returns.

m = nanmean(df);
sigma = nancov(df);
vol = nanstd(df);
cor = corr(df, 'rows', 'pairwise');


save Data/CurveStats m sigma vol cor

Option Impled Volatility

Import a table of option prices for the valuation date and compute implied volatility using Black's model for at-the-money calls and puts.

opt = readtable('Data/NGFuturesOptions.xlsx');

% Only keep at the money puts and calls
opt = filterATMOptions(opt);

% Compute implied volatility
opt.TimeToExp = yearfrac(valDate, opt.Expiration);
opt.IV = blkimpv(opt.UnderlyingPrice, opt.Strike, discRate, opt.TimeToExp,...
                 opt.MidPrice, 1, [], cellstr(opt.Type));

% Average ImpVol (call & put) for each maturity
ivTerm = grpstats(opt,'FutExpiration',@mean,'DataVars','IV');
impVol = ivTerm.mean_IV;

visualizeHistImpVol(vol, ivTerm.mean_IV, ivTerm.FutExpiration(1));

save Data/OptionVol opt ivTerm impVol