The amount of natural gas in storage reflects the ready supply of this commodity, and it is therefore an important factor in pricing natural gas. Often, when estimating a statistical model of natural gas prices, an analyst will estimate a static coefficient for a deseasoned gas storage or weather variable. Doing so implicitly assumes that market participants react identically throughout the year (and over each year) to that variable. Further, it assumes that market participants find the variable no more meaningful in winter, and that they do not adapt to what has occurred in the market. These are unrealistic assumptions of economic behavior. Market participants continuously try to improve their forecasts of market prices, and this likely means they continuously change the scale of their reaction to changes in underlying variables. Moreover, it is likely that natural gas uncertainty is not simply attributable to the regression error term, but also due to changes in how market participants link prices to storage and weather (the coefficients) and the uncertainty in these parameter estimates.

In this analysis, we will model natural gas price returns as a linear function of gas storage and weather variables, and we will allow the coefficients of this function to vary continuously over time. This may be referred to as a time-varying-parameter (TVP) model, and it can be estimated using the Kalman filter. The TVP model allows market participants to adapt their reactions to information contained in variables and affords us an estimate of conditional heteroscedasticity due to both parameter uncertainty and a standard error term.¹¹ 1 Note that an autoregressive conditional heteroscedasticity (ARCH) framework does not specify the source of the conditional heteroscedasticity. We will use this to estimate a time series of the proportion of the total volatility that is attributable to each independent variable, using the parameter uncertainty estimates made available by the Kalman filter. Our results are useful for hedging and derivatives trading in natural gas markets.

Our procedure differs from previous work on volatility spillovers (Diebold and Yilmaz 2009, 2012) in several important ways. First, they consider how total volatility in one variable spills over into another, whereas we measure the effect of volatility induced by uncertainty in the parameters of the function linking the variables. Second, to estimate a time series of volatility spillovers, they must use a rolling window approach. Because our method is based on the Kalman filter, time-varying estimates of parameter uncertainty are built into the model.

Our paper builds on previous analyses of fundamental drivers of volatility in natural gas markets. In particular, Mu (2007) found evidence that weather significantly affects the conditional volatility of natural gas futures returns. The analysis was conducted in a multivariate generalized autoregressive conditional heteroscedasticity (GARCH) framework, which, unlike our Kalman filter approach, does not allow the estimated parameters in the variance equation to vary over time. Geman and Ohana (2009) found a negative correlation between natural gas inventory and price volatility. This negative correlation was particularly pronounced in winter and when inventories were in a period of scarcity. Chiou-Wei et al (2014) showed that deviations in the weekly natural gas storage report from market expectations have a significant effect on natural gas prices. Chan et al (2009) found that both weather and natural gas inventories affect the intensity and size of natural gas futures price jumps. Herbert (1995) showed that trade volume affects the variability of futures volatility, which is evidence of the influence of market microstructure. In an analysis of futures versus spot natural gas prices, Herbert (1993) found that futures prices inaccurately incorporate information into prices. Further, Movassagh and Modjtahedi (2005) found that natural gas futures are biased predictors of future spot prices consistent with normal backwardation. Their results were inconsistent with the bias being explained by the effect of systemic risk on risk premiums in natural gas futures prices.

This paper is organized as follows. In Section 2, we describe the data and provide preliminary evidence for our TVP model. In Section 3, we test for parameter instability and the structural form of time-varying parameters. In Section 4, we introduce our main model and the calculation of conditional variances. In Section 5, we summarize our results and perform diagnostic tests. In Section 6, we highlight the applications of our method. Section 7 concludes.

Our natural gas prices and storage data are from the US Energy Information Administration (EIA) and were gathered from the Administration’s application programming interface using the EIAdata R package (Brigida 2015). Natural gas prices are NYMEX continuous front-month futures prices in USD/MMBtu (EIA API Series ID: NG.RNGC1.W). These contracts are for delivery at the Henry Hub in Louisiana.²² 2 There are many different delivery points throughout the United States, such as Transco Zone 6 (New York City, New York) and Algonquin Citygates (Boston, Massachusetts). Spot prices at these different delivery points are regularly the reference rate for forward and, to a lesser extent, futures contracts. However, the Henry Hub delivery is the most liquid in both the spot and futures markets, and the futures market for delivery at the Henry Hub is where price discovery occurs. We convert these prices into log returns. We convert the storage series into percent deviations from the past five-year norm. The “norm” is defined as the average of the storage value over the same week in each of the previous five years.

Heating degree day ($\mathrm{HDD}$) data is from the US National Oceanic and Atmospheric Administration (NOAA). We use HDD deviations from the norm to measure temperature fluctuations. Consistent with the literature, this deviation is in HDD values and not percent, given the effect of a small denominator in many weeks. This measures both warm () and cold ($\mathrm{HDD}>0$) temperatures.³³ 3 We do not include deviations in cooling degree days because, in general, a large positive deviaton from the norm in $\mathrm{CDD}$ tends to coincide with a large negative deviation from the norm in $\mathrm{HDD}$ (the variables have a less than $-$0.50 correlation coefficient). All data is sampled at a weekly frequency and ranges from January 1, 1999 to June 20, 2014. Summary statistics for our data set are in Table 1. Time series plots of natural gas futures prices and returns are in Figure 1.

We first use a series of ordinary least squares (OLS) regressions to investigate whether the sensitivity of natural gas prices to changes in storage and the weather ($\mathrm{HDD}$) is time varying, and whether the proportion of natural gas volatility attributable to these variables is time varying. To do so, we first calculate the average natural gas log returns ($\mathrm{ng}$), storage deviation ($\mathrm{Stor}$) from a five-year norm, and $\mathrm{HDD}$ deviation for each week of the year. We reject a unit root in each of these variables. We then run fifty-two (one for each week) OLS regressions,

where $w$ denotes a particular week of the year, and $i$ ranges over the fourteen years in our sample. From these regressions, we extract a set of coefficients for each week as well as the proportion of total natural gas volatility attributable to each independent variable.

A plot of the regression coefficients over each week of the year is in Figure 2. Figure 3 shows the proportion of volatility from storage and weather for each week. This offers some evidence that the sensitivity of natural gas returns to storage and the weather varies throughout the year. Further, the proportion of volatility from storage and the weather is generally higher after the fortieth week of the year. These initial results are consistent with our hypothesis, but since the sample is small, they serve to motivate a more thorough examination of these time-varying relationships.

Prior to implementing the time-varying parameter model, we first test for parameter instability. If we can reject stability in the parameters, this is evidence in favor of time-varying parameters. We then investigate the functional form of the state equation.

In this section and onward, we treat our sample as one three-dimensional time series ${({\mathrm{ng}}_t,{\mathrm{Stor}}_t,{\mathrm{HDD}}_t)}^{\prime }$. That is, we no longer divide the sample into fifty-two series, each of which contains the same week over multiple years. Thus, in these parameter stability tests, we are testing if the coefficients on our time series variables are stable from one week to the next.

In order to test for parameter instability, similar to Kim and Nelson (1989), we use a test proposed by Brown et al (1975) to detect departures from constancy in time series regression relationships. The specific test is referred to by the authors as the “homogeneity test” (see Section 2.5. Moving Regressions). The null hypothesis of the test is that the regression parameters are equal at each time point (stable regression coefficients).

The sample period is split into nonoverlapping intervals of arbitrary length $n$, and the “between group over within groups” ratio of mean sum of squares is calculated as the test statistic. Under ${\mathrm{H}}_0$, the test statistic is distributed as $F(kp-k,T-kp)$, where $k$ is the number of regressors, $p$ is the number of intervals and $T$ is the number of observations.

Applying this test to ${\mathrm{ng}}_t={\beta }_0+{\beta }_1{\mathrm{Stor}}_t+{\beta }_2{\mathrm{HDD}}_t+e_t$, for values of $n$ ranging from 20 to 50, we are able to reject our null hypothesis at the 5% level of significance for all $n$ (most $p$-values are below 0.1%). We therefore reject the stability of the regression coefficients.

The preliminary regression results and parameter stability tests are evidence in favor of parameters that vary with time. However, it is likely that market participants (as rational economic agents) use past and present information when deciding on the appropriate present and future coefficients. This motivates a model where the coefficients are allowed to be updated in a Bayesian fashion when new information arrives, much like the views of market participants. An appropriate specification in this case is a time-varying-parameter model, where the parameters are updated using the Kalman filter.

For the main model of our analysis, the measurement equation is

and for coefficient ${\beta }_n$ we let the transition equation take the form of a random walk:

Note ${\mathrm{ng}}_t$ denotes log returns in natural gas futures prices, while $\mathrm{Stor}$ and $\mathrm{HDD}$ represent deviation from normal storage and heating degree days, respectively.

Estimation of the model is done using the Kalman filter and prediction error decomposition. The likelihood function was maximized using the optim function in the R programming language (R Core Team 2014).

This time-varying-parameter model will estimate varying regression coefficients, but it also affords us an estimate of conditional volatility through the conditional variance of forecast errors from the Kalman filter (see Kim and Nelson 1989). The present analysis further decomposes this conditional volatility into the contribution from each factor.

From the Kalman filter, we estimate the conditional variance as

where $x_{t-1}$ is the vector of explanatory variables, $P_{t|t-1}$ is the variance–covariance matrix of the inference on ${\beta }_t$ conditional on information available at time $t-1$ (${\beta }_{t|t-1}$), and ${\sigma }_e^2$ is the variance of the disturbance term.

To calculate the proportion of natural gas volatility attributable to a particular variable (say $\mathrm{Stor}$), we set that variable to zero in $x_{t-1}$ and any row or column in $P_{t|t-1}$ that involves that variable. Denote these as $x_{t-1}^{\mathrm{Stor}}$ and $P_{t|t-1}^{\mathrm{Stor}}$, respectively. We then recalculate

which affords the conditional variance of natural gas prices without that variable. The difference between the full conditional variance and the conditional variance without the variable ($H_{t|t-1}-H_{t|t-1}^{\mathrm{Stor}}$) affords us the conditional variance attributable to that variable.

Below we will highlight two applications of the TVP model of natural gas returns and its accompanying decomposition of natural gas conditional volatility.

In this analysis, we have modeled natural gas returns explicitly, allowing for market participants to learn over time and to react differently to present changes in economic variables. This learning and adaptation, and the attendant parameter uncertainty, constitutes another source of time-varying conditional volatility.

In so doing, we have found evidence of significant variation in the coefficients linking natural gas returns and its underlying fundamental factors. Further, we found that the time series of the Kalman filtered estimates of the $\mathrm{Stor}$ coefficient did not contain a unit root. This implies that we can make inferences about future coefficient values. The modeling and out-of-sample prediction of future $\mathrm{Stor}$ coefficient values would be useful to include in future research. We also found evidence that the weather ($\mathrm{HDD}$) coefficient does contain a unit root and that weather became a less important determinant of natural gas returns in 2007.

In an original application of the TVP model, we decomposed conditional volatility into a time series of each factor contributing to that volatility. This showed that storage is the dominant component of natural gas volatility throughout the year, with weather being the largest contributing factor only during winter periods. Finally, we showed that the results of this analysis have particular applications to hedging and trading in natural gas markets.

The author reports no conflicts of interest. The author alone is responsible for the content and writing of the paper.