Measurement of Connectedness and Frequency Dynamics in Global Natural Gas Markets

We examine spillovers among the North American, European, and Asia–Pacific natural gas markets based on daily data. We use daily natural gas price indexes from 2 February 2009 to 28 February 2019 for the Henry Hub, National Balancing Point, Title Transfer Facility, and Japan Korea Marker. The results of spillover analyses indicate the total connectedness of the return and volatility series to be 22.9% and 32.8%, respectively. In other words, volatility is more highly integrated than returns. The results of the spectral analyses indicate the spillover effect of the return series can largely be explained by short-term factors, while that of the volatility series can be largely explained by long-term factors. The results of the dynamic analyses with moving window samples do not indicate that global gas market liquidity increases with the increasing spillover index. However, the results identify the spillover effect fluctuation caused by demand and supply.


Introduction
The natural gas market has registered an expanding trend, because the fuel transition from coal to natural gas is accelerating in both the industrial and power sectors to reduce greenhouse gas emissions and prevent air pollution, as well as stagnant nuclear power generation. According to BP Statistical Review of World Energy June 2019 [1], the total primary energy consumption in the world grew from 12.4 billion tonnes of oil equivalent (btoe) in 2011 to 13.9 btoe in 2018. Coal consumption was almost constant at 3.6 btoe, while natural gas consumption increased from 2.8 btoe to 3.3 btoe. As a result, the composition ratio of natural gas increased from 22.4% to 23.9%, though the composition ratio of coal decreased from 30.5% to 27.2%. Besides the traditional gas-producing countries in the Middle East and Southeast Asia, others have been increasing their presence as exporting countries, such as Australia (which has been developing large-scale gas wells, including unconventional gas fields) and the United States of America (USA) (which began to export liquefied natural gas (LNG) derived from shale gas). Furthermore, the development of gas fields has been promoted even in Africa recently. Natural gas production skyrocketed from 28 million tonnes of oil equivalent (mtoe), 457 mtoe, and 120 mtoe in 2011 to 112 mtoe, 715 mtoe, and 203 mtoe in 2018, in Australia, the USA, and all of Africa, respectively. On the other hand, consumption is skyrocketing in China and the Middle East, in addition to steady consumption for the members of the Organisation for Economic Co-operation and Development (OECD). Natural gas consumption increased from 24 mtoe, 168 mtoe, and 1142 mtoe in 2011 to 243 mtoe, 476 mtoe, and 1505 mtoe in 2018 in China, all of the Middle East, and the OECD countries.
In the European market, progress in the deregulation of gas business and the decreasing trend of the proportion of long-term contracts linked to crude oil prices has been activating an inter-market the relationship between European and North American markets. Chai et al. [23] study the relationship between the Chinese and the global market. No research has examined the global natural gas market integration, with the exception of that by Silverstovs et al. [24], analyzing the cointegrated relationship between North America, Europe, and Japan based on monthly data from before the shale gas revolution. No studies analyze the spillover effects between North American, European, and Asia-Pacific natural gas markets based on daily data.
We adopt Diebold and Yilmaz's [3] approach to examine spillovers between global natural gas price indexes based on daily data. The correlation coefficient captures only phenomena that do not include the meaning of the relationship. Although the cointegration analysis and the Granger causality test lead to long-term equilibrium and forecast performance, it is difficult to grasp the whole picture at a glance in the case of many variables to be analyzed. Diebold and Yilmaz's [3] approach captures not only pairwise connectedness but also total connectedness. Although Diebold and Yilmaz's [3] approach is not suitable for dynamic descriptions such as multivariate generalized autoregressive conditional heteroscedasticity (GARCH), it is possible to capture dynamic trends with moving window samples. Moreover, we spectrally decompose the Diebold and Yilmaz [3] index by Fourier transform. Baruník and Křehlík [25] utilized the same technique for the first time and later papers applied this approach to the energy market. For instance, Toyoshima and Hamori [26] measure the spillover index for global crude oil markets and decompose it into long-, medium-, and short-term factors. Ji et al. [27] examine the spillover effects between crude oil, heating oil, gasoline, and natural gas in North America and the United Kingdom.
Our contribution to the literature is threefold. The contributions of this paper relate to clarifying the relationships by measuring the connectedness and its frequency dynamics in the global gas market. First, we indicate no progressing integration of the global natural gas market, although we indicate the strong spillover effect between European markets. Moreover, we indicate that the volatility is higher integrated than returns. The Diebold and Yilmaz's [3] approach indicates that the total connectedness of return and volatility is 22.9% and 32.8%, respectively, while each total connectedness is mostly dependent on the pairwise connectedness between European natural gas price indexes. Second, our spectral analyses indicate that long-term factors contribute to volatility spillovers, while short-term factors contribute to return spillovers. We argue that arbitrage might cause short-term return spillovers and long-term memory of volatility might cause long-term volatility spillovers. Finally, our rolling analyses indicate that the above two characteristics continue. Moreover, we argue that regional climate, demand and supply, and incidents might make spillover effects larger or smaller.
The remainder of this paper is organized as follows. Section 2 describes the analyzed data, summary statistics, and preliminary basic analyses. Section 3 explains the adopted methodology. Section 4 presents the empirical results. Section 5 provides a summary of the findings and states conclusions.

Data
We select natural gas price indicators for the North American, European, and Asia-Pacific markets. We use Henry Hub (HH) futures as the North American index. HH is the name of a distribution hub in Louisiana for the natural gas pipeline system in North America, referring to natural gas delivered at that hub. HH futures are listed on the New York Mercantile Exchange. We adopt National Balancing Point (NBP) and Title Transfer Facility (TTF) futures as the European index. NBP futures are futures of natural gas at a virtual trading location in the United Kingdom, being listed on the Intercontinental Exchange (ICE). This market is related to production in the North Sea gas field, which is trading with the European continent by pipelines, and consumption within the United Kingdom. TTF futures are the futures of natural gas at a virtual trading point in the Netherlands, being listed on the ICE. This market reflects the continental European market, which consists of trading by long-haul pipelines and LNG. We utilize the Japan/Korea Marker (JKM) as the Asia-Pacific index, which is associated with the short-term trading market for LNG in the Asia-Pacific region. It is provided by Platts, which assesses benchmark prices in physical energy markets. We obtain the daily data from February 2, 2009 to February 28, 2019 from Bloomberg. Figure 1 shows the time series for these indexes. First, HH fluctuates independently from the other variables. Second, the two European indexes move together. Finally, JKM is linked to the European market over a specific period. As a result of the Great East Japan Earthquake on 11 March 2011, the amount of power generated by nuclear power plants has decreased significantly in Japan. Thus, we can assume that the Asia-Pacific LNG market was tight, and therefore the JKM was soaring. Moreover, these European indexes might have increased as a result of arbitrage trading with the Asia-Pacific LNG market. HH's downtrend around 2012 might be caused by the expectations of increased shale gas production, while the price spike in 2014 was caused by the North American cold wave.
Pacific region. It is provided by Platts, which assesses benchmark prices in physical energy markets. We obtain the daily data from February 2, 2009 to February 28, 2019 from Bloomberg. Figure 1 shows the time series for these indexes. First, HH fluctuates independently from the other variables. Second, the two European indexes move together. Finally, JKM is linked to the European market over a specific period. As a result of the Great East Japan Earthquake on 11 March 2011, the amount of power generated by nuclear power plants has decreased significantly in Japan. Thus, we can assume that the Asia-Pacific LNG market was tight, and therefore the JKM was soaring. Moreover, these European indexes might have increased as a result of arbitrage trading with the Asia-Pacific LNG market. HH's downtrend around 2012 might be caused by the expectations of increased shale gas production, while the price spike in 2014 was caused by the North American cold wave. Table 1 presents the summary statistics of the return and volatility series of each index. There are 2436 observations in each case. The JKM return series has a distribution biased to the right and the other series have distributions biased to the left, as skewness is negative only for the JKM return series. The kurtosis of all series are extremely large. In other words, all series have fat tail distributions. We can reject the hypothesis that each series is normally distributed by the Jarque-Bera statistics calculated from the skewness and kurtosis of each series.    Table 1 presents the summary statistics of the return and volatility series of each index. There are 2436 observations in each case. The JKM return series has a distribution biased to the right and the other series have distributions biased to the left, as skewness is negative only for the JKM return series. The kurtosis of all series are extremely large. In other words, all series have fat tail distributions. We can reject the hypothesis that each series is normally distributed by the Jarque-Bera statistics calculated from the skewness and kurtosis of each series.
As described in Toyoshima and Hamori [26], let x t be a return with mean zero and variance exp(h t ). The SV model can be expressed as follows: where µ, ϕ, and σ τ are the level of log variance, the persistence of log variance, and the volatility of log variance, respectively.

Preliminary Analyses
The condition for the VMA representation of the VAR model is that all variables are stationary. Accordingly, we test for the stationarity status of all series by the augmented Dickey-Fuller (ADF) test. Table 2 presents the results. The ADF test rejects the null hypothesis that all variables have a unit root. Therefore, we can utilize the VMA representation.
We estimate the coefficient of the diagonal Baba, Engle, Kraft, and Kroner (BEKK) GARCH model to calculate the dynamic correlation coefficients. The BEKK model is one of several variations of multivariate GARCH models depending on the formulation of the time-varying variance-covariance matrix. The BEKK model was developed by Engle and Kroner [28], and it guarantees a positive estimated variance. Moreover, the BEKK model can dynamically calculate the correlation coefficient. When there are k variables, the variance-covariance matrix H t of the BEKK model is as follows: where M and N are a k × k matrix, and L is a k × k symmetric matrix. We calculate the correlation coefficients over the entire period, namely from 2 February 2009 to 28 February 2019, to understand the simultaneous relationship between variables as a phenomenon. Tables 3 and 4 present the correlation coefficients for the return and volatility series, respectively. The correlation coefficient between the NBP and TTF is the largest for both the return and volatility series. Furthermore, the correlation coefficients of the volatility series are larger than those of the return series under each combination of variables. This fact indirectly indicates the high sustainability of volatility.  Figure 2 provides the time series of the dynamic correlation coefficients between each pair of return series. Although the NBP and TTF have the highest correlation coefficient in Table 4, even that dynamic correlation coefficient temporarily decreases to 0. Conversely, even if a correlation coefficient over the entire period is low, the dynamic correlation coefficient may temporarily increase.      Figure 3 provides the time series of the dynamic correlation coefficients between each volatility series. We find that the correlation and inverse correlation phases are clear because volatility accumulates at a high rate. However, only the dynamic correlation coefficient between the NBP and TTF continues to be almost +1.

Spillover Index
A correlation coefficient only represents a simultaneous phenomenon between variables. Therefore, our analysis on the relationship between variables needs to go further. We begin by calculating the index proposed by Diebold and Yilmaz [3] which can capture the spillover effect between variables.
We assume the following covariance stationary four-variable VAR ( ): where is the four-dimensional vector of return or volatility, which must be a stationary series; are the 4 × 4 coefficient matrices; is the lag length based on the Schwarz information criterion; and is an independently and identically distributed sequence of four-dimensional random vectors with zero mean and covariance matrix E( ) = Σ.
We can represent the above VAR model in the following VMA: being a 4 × 4 identity matrix with = 0 for 0.
We define the spillover effect from the th to the th market up to -step-ahead as the following equation by using the -step-ahead forecast error variance decompositions: where is the standard deviation of the error term for the th equation and is the selection vector, with one as the th element and zeros otherwise.
We normalize each entry of the variance decomposition matrix by the row sum, that is, 4, as the pairwise connectedness:

Spillover Index
A correlation coefficient only represents a simultaneous phenomenon between variables. Therefore, our analysis on the relationship between variables needs to go further. We begin by calculating the index proposed by Diebold and Yilmaz [3] which can capture the spillover effect between variables.
We assume the following covariance stationary four-variable VAR (p): where x t is the four-dimensional vector of return or volatility, which must be a stationary series; Φ l are the 4 × 4 coefficient matrices; p is the lag length based on the Schwarz information criterion; and ε t is an independently and identically distributed sequence of four-dimensional random vectors with zero mean and covariance matrix E ε t ε t = Σ. We can represent the above VAR model in the following VMA: where A l = p l=1 Φ l A t−l , A 0 being a 4 × 4 identity matrix with A l = 0 for l < 0. We define the spillover effect from the mth to the lth market up to H-step-ahead as the following equation by using the H-step-ahead forecast error variance decompositions: where σ mm is the standard deviation of the error term for the mth equation and e l is the selection vector, with one as the lth element and zeros otherwise. We normalize each entry of the variance decomposition matrix by the row sum, that is, 4, as the pairwise connectedness: We define the sum of pairwise connectedness as total connectedness: The numerator is the sum of the spillover effects, excluding the spillover effects on itself. In other words, the total connectedness means the sum of the relative proportion of the portfolio's response to a shock.
Moreover, we measure the directional spillover effects received by the lth market from all other markets as: Similarly, we measure the directional spillover effects transmitted by the lth market to all other markets as:

Spectral Analysis
Based on the Fourier transform utilized by Baruník and Křehlík [25], we spectrally decompose Diebold and Yilmaz's [3] indexes into short-term (1 to 5 business days) factors, medium-term (6 to 20 business days) factors, and long-term (from 21 business days onward) factors.
The Fourier transform of Equation (2) is as follows: This is the spillover effect from the mth market to the lth market up to H-step-ahead, which is expressed by angular frequency ω.
We define the weighting function, the ratio of ω component to all-frequency components concerning the spillover to the lth variable, as follows: The spillover index from the mth market to the lth market in all bands is expressed as: The spillover index from the mth market to the lth market in band d is expressed as: We convert (θ d ) lm to the relative contribution θ d lm as: We calculate the total spillover index in the d band as: Naturally, this is consistent with Diebold and Yilmaz's [3] spillover index.

Rolling Analysis
It is insufficient to focus on static spillover indicators, which are calculated by the Diebold and Yilmaz [3] and Baruník and Křehlík [25] for the entire period. To capture the dynamics of the spillover effects, we employ the rolling window approach as a sampling method. As shown in Figure 4, this study fixes the moving window sample size to 300 trading days and offsets the window by one business days every time we perform an analysis.
Naturally, this is consistent with Diebold and Yilmaz's [3] spillover index.

Rolling Analysis
It is insufficient to focus on static spillover indicators, which are calculated by the Diebold and Yilmaz [3] and Baruník and Křehlík [25] for the entire period. To capture the dynamics of the spillover effects, we employ the rolling window approach as a sampling method. As shown in Figure 4, this study fixes the moving window sample size to 300 trading days and offsets the window by one business days every time we perform an analysis.

Spillover Index
The spillover analyses results of the return series are reported in Table 5. The spillover index from HH to the other variables and from the other variables to HH are 0.42% and 0.57%, respectively. HH fluctuates independently of the other variables.
The spillover indexes from NBP to TTF and from TTF to NBP are both above 40%. While the total connectedness is 22.9%, the spillover indexes from the other variables to TTF and from TTF to the other variables are 10.7% and 10.8%, respectively. Moreover, the spillover indexes from the other variables to NBP and from NBP to the other variables are 10.5% and 11.1%, respectively. TTF and NBP have a high presence as both the sources and destinations of spillover effects, because these two natural gas markets in Europe are almost integrated, although HH and JKM fluctuate relatively independently. However, only the spillover index from NBP to JKM is not at a negligible level, at 2.65%. Table 5. Spillover index between return series.

Spillover Index
The spillover analyses results of the return series are reported in Table 5. The spillover index from HH to the other variables and from the other variables to HH are 0.42% and 0.57%, respectively. HH fluctuates independently of the other variables.
The spillover indexes from NBP to TTF and from TTF to NBP are both above 40%. While the total connectedness is 22.9%, the spillover indexes from the other variables to TTF and from TTF to the other variables are 10.7% and 10.8%, respectively. Moreover, the spillover indexes from the other variables to NBP and from NBP to the other variables are 10.5% and 11.1%, respectively. TTF and NBP have a high presence as both the sources and destinations of spillover effects, because these two natural gas markets in Europe are almost integrated, although HH and JKM fluctuate relatively independently. However, only the spillover index from NBP to JKM is not at a negligible level, at 2.65%.
The spillover indexes of the volatility series are reported in Table 6. The total connectedness of the volatility series is 32.8%, while the total connectedness of the return series is 22.9%. Risks tend to be transmitted between markets rather than returns. This mutual relationship is similar to the tendency of the return series. However, the characteristic of volatility series is that the influence of HH on the other variables is larger. Additionally, we hardly observe the spillover effects from JKM to the other variables.

Spectral Analysis
The spectral analyses results of the return series spillover are reported in Table 7. The total connectedness from 1 to 5 business days, from 6 to 20 business day, and over 21 business days is 17.3%, 4.12%, and 1.49%, respectively. The short-term factors contribute most to the return spillover. The spillover effects of the return series are mostly explained by events within about one week. Arbitrage might contribute to these short-term spillover effects.  Table 8 presents the spectral analyses results of the volatility series spillovers. The total connectedness from 1 to 5 business days, from 6 to 20 business days, and over 21 business days is 0.10%, 2.37%, and 30.3%, respectively. Contrary to the return series, the long-term factors contribute most to the volatility spillover. Most of the spillover effect of the volatility series is caused by events that occurred more than one month ago. The reason might be the long-term memory of volatility.  Figure 5 shows the total connectedness of the return series and the results of its spectral analyses. The total connectedness has not fluctuated significantly for 10 years and there is no significant change in the composition ratio of the total connectedness. The short-term factors have almost caused return spillovers.   Figure 6 represents the total connectedness of the volatility series and the results of its spectral analyses. As with the static analysis results, the long-term factors almost caused volatility spillovers.

Total Connectedness
Total connectedness spikes in February 2014, February 2018, and November 2018. We assume these spikes are caused by sudden fluctuations in each series. The price of HH spikes in February 2014 because of a great cold wave. The same phenomenon occurs every winter, although the scale is small. Additionally, the JKM, NBP, and TTF prices began to decline at around same time, because the spot market, which has been tight since the Great East Japan Earthquake, shifts and becomes loose. Therefore, the spillover index spikes.
The oversupply strengthens from January to February 2018, because of the expectation to expand USA's LNG export capacity along with the production of crude oil and natural gas. As a result, each gas price index falls sharply. Conversely, in November, the HH increases sharply because the storage stock levels in the USA fall significantly below recent levels. These facts cause the two spikes in 2018.   Figure 7 traces the pairwise connectedness of the return series and the results of its spectral analyses. In all combinations, the spillover effects of the return series continue to depend on shortterm factors. Compared to the total connectedness, which is stable at around 30%, the pairwise connectedness between any variables is low and stable, although the pairwise connectedness from NBP to JKM and from TTF to JKM temporarily reaches around 3%. The mutual spillover effects offset each other because they are almost at the same level. However, when the Asia-Pacific spot market becomes tight or loose, the European market might be relatively stronger than the Asia-Pacific one. Figure 8 shows the pairwise connectedness of the volatility series and the results of its spectral analyses. Under all combinations, the spillover effects of the volatility series continue to depend on long-term factors. The spillover effects from HH to NBP and from TTF to JKM are periodically strong. The HH volatility spikes caused by the cold wave in North America affect the global market. In the European market, the spillover effect from TTF to NBP is stably high. TTF might be more prone to risks than NBP, because the former is more closely linked to intra-European and international trading using pipelines and LNG.  Figure 7 traces the pairwise connectedness of the return series and the results of its spectral analyses. In all combinations, the spillover effects of the return series continue to depend on short-term factors. Compared to the total connectedness, which is stable at around 30%, the pairwise connectedness between any variables is low and stable, although the pairwise connectedness from NBP to JKM and from TTF to JKM temporarily reaches around 3%. The mutual spillover effects offset each other because they are almost at the same level. However, when the Asia-Pacific spot market becomes tight or loose, the European market might be relatively stronger than the Asia-Pacific one.    Figure 8 shows the pairwise connectedness of the volatility series and the results of its spectral analyses. Under all combinations, the spillover effects of the volatility series continue to depend on long-term factors. The spillover effects from HH to NBP and from TTF to JKM are periodically strong. The HH volatility spikes caused by the cold wave in North America affect the global market. In the European market, the spillover effect from TTF to NBP is stably high. TTF might be more prone to risks than NBP, because the former is more closely linked to intra-European and international trading using pipelines and LNG. Figure 7. Pairwise connectedness of return series.

Conclusions
We adopt the approach of Diebold and Yilmaz [3] to examine the spillover effects between global natural gas price indexes. Moreover, we employ the Fourier transform utilized by Baruník and Křehlík [25] for spectrally decomposing Diebold and Yilmaz's [3] index.
We use daily data from 2 February 2009 to 28 February 2019. We employ HH, NBP, TTF, and JKM as the North American, British, Continental European, and Asia-Pacific market price indexes, respectively.
The results of the analyses for the entire period show that the total connectedness of return and volatility is 22.9% and 32.8%, respectively. Volatility is more highly integrated than returns in the global natural gas market. However, total connectedness is mostly dependent on the pairwise

Conclusions
We adopt the approach of Diebold and Yilmaz [3] to examine the spillover effects between global natural gas price indexes. Moreover, we employ the Fourier transform utilized by Baruník and Křehlík [25] for spectrally decomposing Diebold and Yilmaz's [3] index.
We use daily data from 2 February 2009 to 28 February 2019. We employ HH, NBP, TTF, and JKM as the North American, British, Continental European, and Asia-Pacific market price indexes, respectively.
The results of the analyses for the entire period show that the total connectedness of return and volatility is 22.9% and 32.8%, respectively. Volatility is more highly integrated than returns in the global natural gas market. However, total connectedness is mostly dependent on the pairwise connectedness between NBP and TTF. Therefore, we cannot conclude that the global natural gas market is integrated.
The results of the spectral analyses indicate that the spillover effects of the return series almost depend on short-term factors, that is, events within 5 business days. If we further expand the natural gas pipeline network to produce more LNG, the total connectedness of the return series and the ratio of short-term factors should be higher with the increase in arbitrage trading. By contrast, the spillover effect of the volatility series is almost dependent on the long-term factors, that is, accumulated shocks for more than a month. This result represents the long-term memory of volatility.
The results of the dynamic analyses with moving window samples are consistent with the results above. The spillover effects of the return and volatility series tend to be dependent on events within a week and more than a month ago, respectively. However, regional climate, demand and supply, and incidents might make spillover effects unstable. Unfortunately, the spillover index for around 10 years does not show that the liquidity of the global natural gas market tends to increase. For practitioners, this study implies that constantly monitoring market spillovers is significant. The information obtained by this study's approach is also of high value for energy portfolio rebalancing, including derivatives. Moreover, if we monitor the index as a predictive indicator, the information helps with early risk aversion.
Finally, clarification of the spillover effects between the crude oil and natural gas markets will be required for future studies, because many natural gas trading contracts still have their prices linked to crude oil prices.