Abstract
In this work, we examine the relationship between different energy commodity spot prices. To do this, multivariate stochastic models with and without external random interventions describing the price of energy commodities are developed. Random intervention process is described by a continuous jump process. The developed mathematical model is utilized to examine the relationship between energy commodity prices. The time-varying parameters in the stochastic model are estimated using the recently developed parameter identification technique called local lagged adapted generalized method of moment (LLGMM). The LLGMM method provides an iterative scheme for updating statistic coefficients in a system of generalized method of moment/observation equations. The usefulness of the LLGMM approach is illustrated by applying to energy commodity data sets for state and parameter estimation problems. Moreover, the forecasting and confidence interval problems are also investigated (U.S. Patent Pending for the LLGMM method described in this manuscript).
1. Introduction
Understanding the economy evolution in response to structural changes in the energy commodity network system is important to professional economists. The relationship between the different energy sources and their uses provide insights into many important energy issues. The quantitative behavior of energy commodities in which the trend in price of one commodity coincides with the trend in prices of other commodities has always raised the question of whether there is any relationship between prices of energy commodities. If there is any relationship, then what comes to mind is the extent to which one commodity influences the other. Petroleum, natural gas, coal, nuclear fuel, and renewable energy are termed as primary energy components of the energy goods network system because other sources of energy depend on them. Natural gas is usually found near petroleum. Hence, natural gas and crude oil are rivals in production and substitutes in consumption, and energy theory suggests that the two prices should be related. The electric power sector uses primary energy such as coal to generate electricity, which makes electricity a secondary rather than a primary energy source. According to the US Energy Information Administration (EIA), the major energy goods consumed in the United States are petroleum (oil), natural gas, coal, nuclear, and renewable energy. The majority of users are residential and commercial buildings, industry, transportation, and electric power generators. The pattern of fuel usage varies widely by sector [1]. For example, 71% of total petroleum oil provides 93% of the energy used for transportation; 23% of total petroleum oil provides 17% of energy used for residential and commercial use; 5% of total petroleum oil provides 40% of energy used for industrial use; but only 1% of total petroleum oil provides about 1% of the energy used to generate electric power, whereas coal provides 46% of the energy used to generate electric power and natural gas provides 20% of the energy used to generate electric power. This analysis suggests that the strength of interactions between coal and electricity will be stronger than that of crude oil and electricity, or natural gas and electricity.
Energy price forecasts are highly uncertain. We might expect the price of the natural gas and crude oil to follow the same trend because they are often found mixed with oil in oil wells and also of the fact that natural gas is often used in petroleum refining and exploration. Recently, Ramberg and Parsons [2] showed that the cointegration relationship between natural gas and crude oil does not appear to be stable through time. They claimed that though there is cointegration between the two energy prices, there are also recurrent exogeneous factors such as seasonality, episodic heat waves, cold waves, and supply interruption from the hurricane affecting the trends in the prices. Brown and Yücel [3] also observed that the price of natural gas pulled away from oil prices in 2000, 2002, 2003, and late 2005. Oil prices are influenced by several factors, including some that have mainly short-term impacts and other factors, such as expectations about the future world demand for petroleum, other liquids, and production decisions of the Organization of the Petroleum Exporting Countries (OPEC) [1]. Supply and demand in the world oil market are balanced through responses to price movement with considerable complexity in the evolution of underlying supply and demand expectation process. For petroleum and other liquids, the key determinants of long-term supply and prices can be summarized in four broad categories [1]: the economics of non-OPEC supply, OPEC investment and production decisions, the economics of other liquids supply, and world demand for petroleum and other liquids.
An understanding of how changes in the price of one energy commodity are expressed in terms of other energy commodity is needed. This would prove to be useful in predicting price behavior over the long run and further facilitates profit-maximizing processes. To check if there is really indeed a relationship between energy commodities, the need to be able to create a model which explains such commodity price relationship over short- and long-time intervals arises. The relationships between energy commodities have been addressed in [2–15]. The error correction model [2–5, 7, 8] is the most commonly used model by authors to describe the relationship between energy commodities. Moreover, Hartley et al. [7] have concluded that the U. S. natural gas and crude oil remain linked in their long-term movements. In addition, it is exhibited that there is strong evidence of a stable relationship between these two energy commodities which are affected by short-run seasonal fluctuations and other factors. The rule of thumb [7] has long been used in the energy industry to relate the natural gas prices to crude oil prices. The rule denoted by the 10-to-1 rule states that the price of natural gas is one-tenth of the price of crude oil prices. Similarly, the 6-to-1 rule states that the price of natural gas is one-sixth of the price of crude oil. It has been examined by Brown and Yücel [3] that these two rules do not perform well when used to assess the relationship between US natural gas price and West Texas Intermediate (WTI) crude oil price for the past years. Moreover, their error-correcting model specifies the relationship between the two commodities. Using this model, they show that when certain factors are taken into account, movements in crude oil prices can shape the price of natural gas. Vezzoli [9] in his work applies an ordinary least squares (OLS) regression on the log of natural gas and crude oil prices. Using this model, he was able to conclude that there is a relationship between natural gas and crude oil prices. Bachmeir and Griffin [4] showed that crude oil, coal, and natural gas in the United States have weak cross-cointegration using the error correction model. Ramberg and Parsons [2] show that any simple formula between natural gas and crude oil prices will leave a portion of the natural gas price unexplained. Furthermore, the relationship between natural gas and crude oil using a vector error correction model [2, 3] under the cointegration between the two energy commodities and other factors such as recurrent exogeneous factors are presented. Villar and Joutz [10] list some economic factors linking natural gas and crude oil prices, while testing for market integration in the United Kingdom during the time when natural gas was deregulated. Asche et al. [6] have integrated the prices of the energy commodities: natural gas, electricity, and crude oil.
Most of this work is centered around the relationship between prices of energy commodities. We are interested in an interdependence of certain energy commodities. Moreover, a nonlinear multivariate interconnected stochastic model of energy commodities and sources with and without external random intervention processes is developed. Also, parameter and state estimation problems of continuous time nonlinear stochastic dynamic process are motivated to initiate an alternative innovative approach. This led to the introduction of the concept of statistic processes, namely, local sample mean and sample variance which further led to the development of an interconnected discrete-time dynamic system of local statistic processes and its mathematical model discussed in Otunuga et al. [16, 17]. This paved the way for extending the LLGMM [16] to a multivariate local lagged adapted generalized method of moments case. The parameters in the multivariate model are estimated using the LLGMM method. These estimated parameters help in analyzing the short- and long-term relationships between the energy commodities. It has been shown in Appendices A, B, C, and D of the work of Otunuga et al. [16] that the LLGMM parameter estimation scheme performs better than existing nonparametric statistical methods. A performance comparison (with appropriate statistical indices) of the LLGMM method with existing orthogonality condition based on generalized method of moments (OCBGMM) and aggregated generalized method of moments (AGMM) methods using the energy commodity stochastic model is discussed in their work. The method is applied to estimate time-varying parameter estimates in a stochastic differential equation governing energy commodities, stock price processes, and transmission of infectious diseases in the work of Otunuga et al. [16], Otunuga [18], and Mummert and Otunuga [19], respectively. Algorithm for implementing the LLGMM approach is presented in Otunuga et al. [17]. In this work, the usefulness of this approach is further illustrated by applying the technique to study the relationship between three energy commodity data sets: Henry Hub natural gas, crude oil, and coal data sets for state and parameter estimation problems. Interested readers are directed to the work of Otunuga et al. [16–18, 20] to read more about the LLGMM parameter estimation procedure. Moreover, the forecasting and confidence interval problems are also investigated.
The organization of this paper is as follows.
In Section 2, we derive a multivariate stochastic dynamical model for energy commodities. In Section 3, the multivariate model derived in Section 2 is validated. Due to random intervention, we introduce interventions described by a continuous jump process in our model in Section 4. In Section 5, we discuss the discrete-time dynamic model for sample mean and covariance processes with jump using the work of Otunuga et al. in [16, 17]. In Section 6, we discuss about parametric estimation techniques. We construct an observation system from a nonlinear stochastic functional differential equation. In Section 7, using the method of moments [21], in the context of lagged adaptive expectation process [22], we briefly outline a procedure to estimate the state parameters for the dynamical model with jump and the model without jump locally. Moreover, the usefulness of computational algorithm is illustrated by applying the procedure to test for the relationship between Henry Hub natural gas, crude oil, and coal for the state and parameter estimation problems. In Section 8, the forecasting and confidence interval problems are also addressed.
2. Model Derivation
Let be a vector of energy commodity prices considered to have long-run or short-run relationship with each other, with being the price of the th energy commodity at time . The economic principles of demand and supply processes suggest that the price of an energy commodity will remain within a given finite expected lower and upper bounds. We define and as the expected upper and lower limits of the th energy commodity spot prices, respectively. In the absence of interactions between the energy commodities , , where the market potential for the th commodity per unit of time at time can be characterized by . This simple idea leads to the following economic principle regarding the dynamic of the price of the th energy commodity. The change in spot price of the th energy commodity over the interval of length is directly proportional to the market potential price, that is,
This implies that for some constant . From this deterministic mathematical model, if , we note that as the price falls below the expected price , the price of the th commodity rises, and as the price lies above , there is a tendency for the price to fall. Similar argument follows if . Hence is the equilibrium state of (3).
In a real world situation, the expected upper price limit of the th commodity is not a constant parameter. It varies with time, and moreover it is subjected to random environmental perturbations. Therefore, we consider where is a white noise process that characterizes the measure of random fluctuation of the upper price limit of the th commodity; here stands for the mean of the energy spot price process of the th commodity at time . It is further assumed that is governed by a similar dynamic forces described in (3), that is, where is defined as the mean reversion rate of the mean of the th commodity. By following the argument used in (4), we incorporate the effects of random environmental perturbations into the lower limit of the mean of the th commodity: where , and is a white noise process describing the measure of random influence on the mean price of the th commodity.
Substituting expressions in (4) and (6) into (3) and (5), respectively, we obtain
In the absence of interactions and using (7), the system of stochastic model for isolated expected spot and spot prices processes are described by the following nonlinear system of stochastic differential equations: where and are nonnegative for .
In the presence of interactions, for each , we consider both deterministic and stochastic interaction functions. For each , we define the th aggregate interaction functions and for the th mean energy spot price process and the energy spot price process in energy commodity market network system, respectively. Moreover, we assume that these functions have the following structural forms: where and are elements of interconnection matrices denoted by and , respectively. In (10), and represent a degree of interaction of the th commodity with the th commodity in the energy commodity market network system.
For the matrix , with fixed if the th commodity in the energy market network system does not influence the th commodity. Likewise, for the matrix , with fixed , if the th commodity in the energy market network system subcomponent of is totally unaffected by the influence of the th commodity. Finally, we introduce interactions in the diffusion coefficients with respect to the th commodity of the energy market network system under random environmental perturbations as and for each . The diffusion part is of the form where and are -dimensional white noise processes; stands for dot product.
We assume that the interaction functions (10) and (11) have the following forms: where , are defined in (10), , and , and ; , and are column vectors; and are block diagonal matrices; , . We also assume that , , , and satisfy the local Lipschitz condition. This assumption implies that , , , and also satisfy local Lipschitz condition.
Thus, the interconnected energy commodity network system is described by where the parameters ; ; ; ; ; ; ; and for , ; ; ; for , and are -dimensional independent Wiener processes defined on a filtered probability space (); for , , and for , ; the filtration function is right continuous; for each , each contains all -null sets in ; the -dimensional random vectors and are measurable.
The network system of stochastic differential equations in (13) can be written as follows: where and
and are block matrices; , ; and , are a block matrix with each entries having order .
In the next section, we outline the model validation problems of (14), namely, the existence and uniqueness of solution process.
3. Mathematical Model Validation
Here, we validate the mathematical model derived in Section 2. We note that the classical existence and uniqueness theorem [23–25] is not directly applicable to (14). We need to modify the existence and uniqueness results. The modification is based on Theorem 3.4 [23]. We show the global existence of the solution process of a system of differential equations (14). We note that the rate functions , , , and in stochastic system of differential equations (14) do not satisfy the classical existence and uniqueness conditions [23]. However, these rate functions do satisfy the local Lipschitz condition. Therefore, we construct sequences of functions for the drift and diffusion coefficients of interconnected dynamic system (14) so that the classical conditions for existence and uniqueness theorem are applicable. The construction of modification scheme is as follows: First, we define a cylindrical subset of for , , where is an - dimensional sphere with radius defined by for any . We note that is inscribed in -dimensional parallelepiped in .
The developed stochastic network model (14) can be written as: where
Here, for each and , we define . Hence, .
Remark 1. We observe that satisfies the global Lipschitz condition on with a Lipschitz constant 1. Using this, together with the local Lipschitz condition assumption on the drift and diffusion coefficients of stochastic differential equations (14), the modified rate coefficient functions in (18) satisfy the classical existence and uniqueness conditions [26, 27]. Thus, its solution is denoted by (), for . Moreover, () is nonnegative whenever , .
Now, we apply Theorems 3.4 and 3.5 of [26] in the context of the modified system of stochastic differential equations (18) and Remark 1 to establish the global existence of solution of stochastic differential equations in (13). For this purpose, we outline the argument used in the proof of these theorems. In addition to the local Lipschitz conditions on the drift and diffusion coefficients, we further impose the following hypothesis on the coefficients:
(H1) where for , , are nonnegative; , , . From (18), we further remark that the dynamic of the mean spot price vector is decoupled with the dynamic of spot price . Now, we first apply Theorems 3.4 and 3.5 of [23] in the context of modified system of stochastic differential equations (18) and hypothesis (H1) to establish the global existence of solution of the completely decoupled subsystem of stochastic differential equations in (18). For this purpose, we outline the argument used in the proof of these theorems.
Definition 2. Let be the first exit time of the solution process from the set . Define to be the (finite or infinite) limit of the monotone increasing sequence as .
We wish to show that
In the following, we present a result that is parallel to Theorem 3.5 [26] in the context of the completely decoupled subsystem of stochastic differential equation (14). For this purpose, it is enough to exhibit the global existence result for the transformed system (18).
Lemma 3. For , and , let be the solution of the completely decoupled subsystem of (18), and let the hypothesis (H1) be satisfied. Let be a function defined on into defined by Then, there exists some constant such that where is the differential operator with respect to (14); .
Moreover, the global existence of solution of the completely decoupled subsystem of (14) is follows.
Proof. It is obvious that on . In fact, , , and exist and are continuous functions defined on . Moreover, the operator with respect to the completely decoupled component is as follows:
where
Furthermore, as .
To show that , we define a function
It is obvious that . By defining ; for ; and imitating the argument of Lemma 4 [23], we have
Hence,
The global existence and uniqueness of the solution of the first component of (18) follows by letting . Hence, from this and the fact that the solution process of transformed system (18) is almost surely identical with the solution process of the original system (14), we conclude the global existence and uniqueness of (14).
Following the idea of Lemma 3, we present a global existence and uniqueness of solution of the system of stochastic differential equations governing the subsystem in (14).
Lemma 4. For , and , , let be the solution of the system of stochastic differential equations governing the subsystem described in (18), and let the hypothesis (H1) be satisfied. Let be a nonnegative function on into defined by Then, there exist a constant such that where is the differential operator with respect to (14); .
Moreover, the global existence of the solution of the system of stochastic differential equations governing the subsystem described in (14) is as follows.
Proof. It is obvious that on . In fact, , , , and exist and are continuous functions defined on . Moreover, the operator with respect to the system of stochastic differential equations governing the subsystem described in (14) is obtained as follows: where Furthermore, as . Following the final argument used in proving the global existence of in Lemma 3, we conclude that there exists a unique global solution to the interconnected system of stochastic differential (14).
In the next section, we discuss the case where we incorporate jump process into the system ().
4. Energy Commodity Model with and without Jumps
Due to random interventions that affects the price of energy commodities, we introduce random perturbations described by a continuous jump in our model. We follow the approach discussed in [28, 29] where a class of stochastic hybrid dynamic processes is investigated.
Let be the number of jumps on the interval , for . For , let be the jump times over a time interval such that , where denotes the time at which the th jump occurred in the system . For , no jump has occurred on the interval . We denote the th subinterval by . Knowing the global existence and uniqueness solution process of the system (14) on the interval , in Section 3, for , where if , and if and , we consider the solution process on each subinterval , between jumps. For and , we have and . In this case, we consider the solution process on . The interconnected system is governed by the following system of continuous time stochastic process: where is the solution of the system of equation (34); for and , and are jump coefficient matrices. These jump coefficients are estimated by ; for . Following the approach in [28, 29], the solution of (34) takes the following representation:
Remark 5. In case of no jump, that is , (34) and (36) reduce to respectively.
5. Discrete-Time Dynamic Model for Local Sample Mean and Covariance Processes
In this section, we extend the discrete-time dynamic model for local sample mean and variance processes given in Otunuga et al. [16, 17] to a multivariate case including the jump process. The development of idea and model of statistic for mean and covariance processes was motivated by the state and parameter estimation problems of the continuous time nonlinear stochastic dynamic model of the energy commodity market network. For this purpose, we need to introduce a few definitions and notations. For more information, we direct the readers to the work of Otunuga et al. [16, 17].
For each , let and be finite constant time delays such that . If , then there is no jump and we simply write and . Here, characterizes the influence of the past performance history of state of dynamic process. describes the reaction or response time delays. In general, these time delays are unknown and random variables. These types of delays play an important role in developing mathematical models of continuous time [30] and discrete-time [22, 31] dynamic processes. Based on the practical nature of data collection process, it is essential to either transform these time delays into positive integers or design a suitable data collection schedule or discretization process. For this purpose, for each , we describe the discrete version of time delays of and as follows:
Here, is the integer part of a real number, and for the sake of simplicity, we assume that .
Definition 6. Let be a continuous time multivariate stochastic process defined on an interval into , for some . For , let be an increasing subsigma algebra of a complete probability space for which is measurable. For each , let and be a partition of and , respectively. The partition is derived by decomposing the partition . For each , the partitions and are defined as follows: where represents time in the interval , , and .
Remark 7. We define with . For , we note that we can write as . From this, it follows directly that the jump times are contained in . For any , , we have . Hence, there is a relationship between elements of with that is described by , for . In fact, the relationship between a set of jump times and the partition defined in (40) is as , where is the size of partition of the subinterval . It follows that . Using these facts, and noting that if , then , , , , , , and (40) reduces to For each , let be a finite sequence corresponding to the stochastic process and partition defined in (41). We simply write . We further recall that is measurable for . We also recall the definition of forward time shift operator [16, 17, 32]:
Definition 8. For and , each and each , a partition of closed interval is called local at time , and it is defined by Moreover, is referred to as the point subpartition of the partition in (41) of the closed subinterval of .
Remark 9. We note that for , that is, no jump, we have , , , and , where is referred as the point subpartition of partition in (40) of the closed subinterval of for .
Definition 10. For each , , and , a local finite sequence at of size is the restriction [33] of to in (43). This restriction sequence is defined by
As varies from 2 to , the corresponding respective local sequence at varies from to . As a result of this, the sequence defined in (44) is also called a local moving sequence. Furthermore, the average corresponding to the local sequence in (44) is defined by
The average/mean defined in (45) is also called the local average/mean. For and , the local covariance matrix corresponding to the local sequence in (44) is defined by where , is the local sample covariance statistic between and at described by
Denoting for and , the discrete-time interconnected dynamic model of local sample mean and variance processes (DTIDMLSMVSP) [16, 17] is given by where is the order of the system and
The proof of this is given in [16].
6. Parametric Estimation
In this section, we consider a parameter estimation problem in drift and diffusion coefficients of (34). This is achieved by utilizing the lagged adaptive process [22] and the interconnected discrete-time dynamics of local sample mean and variances statistic processes model (48). For each , we consider a general interconnected hybrid system described by the system of stochastic differential equations: where is the jump coefficient matrix; the jump times ’s are defined in (34). For each , the estimate of the jump coefficient is given by .
Let , and its partial derivatives , , and exist and are continuous on each interval . We apply Itô-Doob stochastic differential formula [34] to , and we obtain where the operator is defined by
For (50) and (51), we present the Euler-type discretization scheme [23]:
Define as the filtration process up to time . With regard to the continuous time dynamic system (50) and its transformed system (51), the more general moments of are as follows: where , and stands for the transpose of the matrix. From (53) and (54), we have
This provides the basis for the development of the concept of lagged adaptive expectation process [22] with respect to continuous time stochastic dynamic system (51). This indeed leads to a formulation of local generalized method of moments at .
In the following, we state a result that exhibits the existence of the solution of the system of nonlinear algebraic equations. For the sake of easy reference, we shall state the implicit function theorem without proof.
Theorem 11. Implicit Function Theorem [33]. Let be a vector-valued function defined on an open set with values in . Suppose on . Let () be a point in for which and for which the determinant of the Jacobian matrix det . Then, there exists a -dimensional open set containing and unique vector-valued function , defined on and having values in , such that on , , and for every .
6.1. Special Case: Illustration
For each and each , we consider a special case of (14).
Here, , , , , , and are the system parameters on the jump subinterval ; , , , , and are positive; and for , ; , are nonnegative. and are independent standard Wiener process on a filtered probability space () with the properties described in (14). It follows that the interconnected system of stochastic differential equations (56) has parameters. Also, and likewise,
Remark 12. For the case , (56) is reduced to
where for , the parameters , , , , , and are the system parameters on the interval ; , , , , and are positive; and for , ; , are nonnegative. For each , we pick a Lyapunov function
in (51) for (56). Using the Itô differential formula [24], we have
By setting , , and , the combined Euler discretized scheme for (61) is
where , are given finite sequence of measurable random vectors and are independent of , , respectively. We define and .
Applying conditional expectation to (62) with respect to , we obtain
where is the filtration up to time . From (62), (63) and (64) are reduced to
Equation (65) provides the basis for the development of the concept of lagged adaptive expectation process [22] with respect to continuous time stochastic dynamic systems (56) and (61).
For , applying the lagged adaptive expectation process [22] and Definitions 6, 8, and 10 and using (63), (64), and (65), we formulate a local observation/measurement process at as an algebraic functions of local functions of restriction of the finite sample sequence and to subpartition in Definition 8:
For each and each , we define by
For every , we have
We define , , , and . Thus, provided that the determinant of each of the Jacobian matrices , , and are not zero, by the application of Theorem 11 (Implicit Function Theorem), we conclude that for every nonconstant -local sequence and , there exist a unique solution (, , , , , ) of system of algebraic equations (70) as a point estimates of , , , , , , respectively. We illustrate this approach using energy commodities natural gas, crude oil and coal [35–37] in the next section.
6.2. Illustration: Application to Energy Commodity
In this subsection, we present an illustration of the derived model on the natural gas, crude oil and coal [35–37]. Here, and . Moreover, (56) reduces to
For each , following the argument used in Section 6.1, we have
We also have , , , and . For each , it follows that the determinant of the Jacobian of , that is,
, is not zero provided that all parameters are not zero or provided the sequence is neither zero nor constant. Likewise, determinants of the Jacobians , and are nonzero if , , and are not zero or provided the sequence and are neither zero nor constant for .
Remark 13. If the sample is a constant sequence, it follows from (62) (q = 1) and the fact that and that we can set . It also follows from (66) that .
7. Computational Algorithm
In this section, we outline computational, data organizational, and simulation schemes. We introduce the ideas of iterative data process and data simulation time schedules in relation to the real-time data observation/collection schedule. For the computational estimation of continuous time stochastic dynamic system state and parameters, it is essential to identify an admissible set of local conditional sample average and sample covariance parameters, namely, the size of local conditional sample in the context of a partition of time interval . Moreover, the discrete-time dynamic model of conditional sample mean and sample covariance statistic processes in Section 5 and the theoretical parameter estimation scheme in Section 6 motivates to outline a computational scheme in a systematic and coherent manner. A brief conceptual computational scheme and simulation process summary is described below.
7.1. Coordination of Data Observation, Iterative Process, and Simulation Schedules
Without loss of generality, we assume that the real data observation/collection partition schedules , are defined in (41). We present definitions of iterative process and simulation time schedule following the definitions in Otunuga et al. [16].
Definition 14. The iterative process time schedule in relation with the real data collection schedule is defined by
where is a forward shift operator [32].
The simulation time is based on the order of the time series model of local conditional sample mean and covariance processes in (48).
Remark 15. For the case where , we have , where is defined in (40). This is the iterative time schedule in the absence of jumps.
Definition 16. The simulation process time schedule in relation with the real data observation schedule is defined by
Remark 17. For each , the initial times of iterative and simulation processes are equal to the real data times and , whenever and , respectively. The iterative process and simulation process times with jump are and , , respectively.
7.2. Conceptual Computational Parameter Estimation Scheme
For the conceptual computational dynamic system parameter estimation, we need to introduce a few concepts of local admissible sample/data observation size local admissible conditional finite sequence at and local finite sequence of parameter estimates at .
Definition 18. For each , and , we define local admissible sample/data observation size at as , where Moreover, is referred as the local admissible set of lagged sample/data observation size at .
Remark 19. We note that if , , the point . Thus, (76) is reduced to
Definition 20. For each , in Definition 18 and , a local admissible lagged adapted finite restriction sequence of conditional sample/data observation at time to subpartition of in Definition 8 is defined by . Moreover, a class of admissible lagged adapted finite sequences of conditional sample/data observation of size at is defined by
In the case of energy commodity model, for each , , we find corresponding local admissible adapted finite sequence of conditional sample/data observation at , . For , using this sequence and solutions of (70), we compute for . This leads to a local finite sequence of parameter estimates at defined on as follows:
The above defined collection is denoted by
7.3. Conceptual Computation of State Simulation Scheme
For the development of a conceptual computational scheme, we need to employ the method of induction. The presented simulation scheme is based on the idea of lagged adaptive expectation process [17]. For , an autocorrelation function (ACF) analysis [21, 32] performed on suggests that the interconnected discrete-time dynamic model of local conditional sample mean and sample variance statistics in (48) is of order . In view of this, we need to identify the initial data. We begin with a given initial data , ,, as defined in (45), (46), (47), and (48).
Let be a simulated value of at time corresponding to an admissible sequence . For , and , the simulated value is generated from the discretized Euler scheme (62) as follows:
To find the simulated value and , we need to estimate and by first simulating
and as follows:
From this, we calculate and as:
Thus, and . Let (, ) be a - local sequence of simulated values corresponding to -admissible lagged adapted finite sequence of conditional observation belonging to , and corresponding term of sequence.
(). Thus, for each , (, ) are the finite sequence correspondence of simulated values of at .
7.4. Mean-Square Suboptimal Procedure
To find the best estimate of using a local admissible finite sequence , we need to compute a finite sequence of quadratic mean square error corresponding to . The quadratic mean square error is defined below.
Definition 21. For each , the quadratic mean square error of relative to each member of the term of local admissible sequence of simulated values is defined by
For any arbitrary small positive number and for each time , to find the the best estimate from the admissible simulated values of simulated sequence of for , we determine the following sub-optimal admissible set of -size local conditional sample
Among these collected values, the value that gives the minimum for are recorded as . If more than one value exist, then the largest of such ’s is recorded as . If condition (86) is not met at time , the value of where the minimum is attained is recorded as . The level sub-optimal estimates of the parameters are recorded as . Finally, the simulated value , at time with is now recorded as the best estimate for and . The value , is called the sub-optimal simulated value of and of and at .
7.5. Illustration: Application of Conceptual Computational Algorithm to Energy Commodity Data Set
In this subsection, we apply the above conceptual computational algorithm to study the relationship between three energy commodities by setting in (56). The three energy commodities are daily Henry Hub Natural gas data set, daily crude oil data set, and daily coal data set for the period of 05/04/2009-01/03/2014, [35–37]. Thus, for each pair , , and , the drift and diffusion coefficient function of the stochastic dynamic equation governing , for has 4 and 3 parameters each to be estimated, respectively. Thus, there are 42 parameters to be estimated in total. Using , , for each , the level sub-optimal estimates of parameters , , , , , , , at each real data times are exhibited below.
7.5.1. Illustration: Relationship between Natural Gas, Crude Oil, and Coal: Without Incorporating Jump Process
In this subsubsection, we analyze the relationship between natural gas, crude oil, and coal without the jump process. For , the stochastic dynamic system governing the three energy commodities is described in (59) of Remark 12. Here, () denotes the mean spot and the spot price process of natural gas, () denotes the mean spot and the spot price process of crude oil, and () denotes the mean spot and the spot price process of coal.
Using the discretized scheme (82), we apply the above conceptual computational algorithm for the real-time data sets namely daily Henry Hub natural gas data set, daily crude oil data set, and daily coal data set for the period of 05/04/2009-01/03/2014, [35–37]. Using , , , and , the -level suboptimal estimates of the parameters at each real data times are described below for each commodity data sets.
The parameters corresponding to the natural gas data set are , , , , , , , , , , , , , and . The parameters corresponding to the crude oil data set are , , , , , , , , , , , , , and . The parameters corresponding to coal data set are , , , , , , , , , , , , , and .
The following table gives the drift coefficient’s parameter estimates , , , , , , , , , , , and for the decoupled dynamical system for in the case where jump is not incorporated into the dynamical system.
Table 1 shows the estimates of the suboptimal size , , the parameters , , , , , , , , , , , and for each of the energy commodity data sets. Moreover, and the initial real data time is .
Figure 1 shows the drift coefficient’s parameter estimates , , and for the decoupled dynamical system for in the case where jump is not incorporated into the dynamical system.
(a)
(b)
(c)
Figures 1(a)–1(c) are the graphs of , , and against time for the daily Henry Hub natural gas price [36], daily crude oil price [37], and daily coal price [35] data set, respectively. By plotting the real data sets, it is easily seen that the graphs of , , and are similar to the graph of the Henry Hub natural gas, crude oil, and coal data set, respectively. We expect this to happen because , , are the equilibrium spot price processes of (5). The sample means of the real spot price data sets for natural gas, crude oil, and coal are given by 3.7218, 88.5620, and 16.5543, respectively. It can be seen from Figures 1(a)–1(c) that the graphs of , , and move around the mean values 3.7218, 88.5620, and 16.5543, respectively. This analysis shows that the parameters , , are statistic processes for the respective mean of the data sets at time .
The graph of the interaction parameters , , , , , , , , and for the decoupled dynamical system for (with no jump incorporated into the dynamical system) are given in Figure 2 below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Figures 2(a)–2(i) show the graph of the suboptimal interaction coefficient parameters , , , , , , , , and . The interaction coefficients , are negligible, because each estimate is . Thus, this shows that the model describing the mean spot price, , is mainly characterized by the market potential , .
Table 2 below shows the estimates of the diffusion coefficient’s parameters , , , ,, , , , and for the decoupled dynamical system for in the case where jump is not incorporated into the dynamical system.
The graph of the diffusion coefficient’s parameter for the decoupled dynamical system for without jump incorporated into the dynamical system are given in Figure 3 below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Figures 3(a)–3(i) show the graphs of the suboptimal interaction measure of fluctuation coefficient parameters , , , , , , , , and .
Table 3 gives the drift coefficient’s parameter estimates , , , , , , , , , , , and for the dynamical system for without jump.
Table 3 shows the estimates of the suboptimal size , and the parameters , , , , , , , , , , , and , for each of the energy commodity data sets. Moreover, , and the initial real data time is .
In the following, we give the graph of the drift coefficient’s parameters with estimates in Table 3 in Figure 4 below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Figures 4(a)–4(i) show the graphs of the suboptimal interaction coefficient parameters , , , , , , , , and . According to (58), the estimate , , is positive if commodity is cooperating with commodity , and negative if commodity is competing with commodity . There is no interaction between the two commodities if . It is apparent from the graph of that coal is competing with natural gas because the estimates of are mostly negative. Also, it is apparent that natural gas and crude oil are either cooperating or competing, depending on the time period.
Figures 5(a)–5(c) are the graphs of , , and against time for the daily Henry Hub natural gas price data set [36], daily crude oil price data set [37], and daily coal price data set, respectively.
(a)
(b)
(c)
Table 4 gives the suboptimal estimates of the parameters , , , , , , , , and for each of the energy commodity data sets.
Figures 6(a)–6(i) are the graphs of , , , , , , , , and against time for the daily Henry Hub natural gas price data set [36], daily crude oil price data set [37], and daily coal price data set, respectively.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Table 5 shows the real and simulated estimates for the spot price processes , .
Figure 7 below shows the graph of the real and simulated prices for natural gas, crude oil, and coal data set.
(a)
(b)
(c)
Figures 7(a)–7(c) show the graph of the real and simulated spot prices for the daily Henry Hub natural gas data set [36], daily crude oil data set [37], and daily coal data set [35], respectively. The red line represents the real data set , while the blue line represent the simulated data set . Here, we begin by using a starting delay of . The simulation starts from . It is clear that the graph fits well. To reduce magnitude of error, we increase the magnitude of time delay. We later compare this result with the case where jump is incorporated into the system.
7.5.2. Relationship between Natural Gas, Crude Oil, and Coal: With Jump Incorporated
In this subsubsection, we analyze the relationship between natural gas, crude oil, and coal with the jump process. Here, we apply the above conceptual computational algorithm described in this section for the real-time data sets, namely, daily Henry Hub natural gas data set, daily crude oil data set, and daily coal data set for the period of 05/04/2009-01/03/2014, [35–37]. For , , we use , , , and . The -level suboptimal estimates of the parameters at each real data times are described below for each commodity data sets.
The parameters corresponding to the model governing natural gas price data set are , , , , , , , , , , , , , and . The parameters corresponding to the model governing crude oil price data set are , , , , , , , , , , , , , and , while the parameters corresponding to the model governing coal price data set are , , , , , , , , , , , , , and .
For the sake of simplicity and in order to be able to compare our results in this subsection with the results in Subsection 7.5.1, for each , we rewrite the parameters , , , , , and after they have been estimated as , , , , , and .
We give some results for the jump times of the system .
Table 6 shows some result for the jump times of the system (). These results are derived by recording the times at which an entry of a commodity differ by twice the standard deviation or more from the mean of that commodity. These times are now combined into a single array and sorted out in an increasing order.
We give the estimate of the jump coefficient matrices and defined in (34) in Table 7 below.
Table 8 gives the estimates , , , , , , , , , , , and for the decoupled dynamical system for with jump.
Figure 8 shows the parameter estimates , , and for the decoupled dynamical system for with jump.
(a)
(b)
(c)
Figures 8(a)–8(c) are the graphs of , , and against time for the daily Henry Hub natural gas price data set [36], daily crude oil price data set [37], and daily coal price data set, respectively. By plotting the real data sets, it is easily seen that the graphs of , , and are similar to the graph of the Henry Hub natural gas, crude oil, and coal data set, respectively. We expect this to happen because , , are the equilibrium spot price processes of (5).
The graph of the interaction parameters , , , , , , , , and for the decoupled dynamical system for with jump and estimates in Table 8 are given in Figure 9 below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Figures 9(a)–9(i)) show the graph of the suboptimal interaction coefficient parameters , , , , , , , , and . The interaction coefficients , , are negligible, because each estimate is . Thus, this shows that the model describing the mean spot price, , is mainly characterized by the market potential , , .
Table 9 below shows the estimates of the diffussion coefficient’s parameters , , , ,, , , , and for the decoupled dynamical system for with jump.
The graph of the diffusion coefficient’s parameter for the decoupled dynamical system for with jump are given in Figure 10 below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Figures 10(a)–10(i) show the graph of the - sub-optimal interaction measure of fluctuation coefficient parameters , , , , , , , , and .
Table 10 gives the drift coefficient’s parameter estimates , , , , , , , , , , , and for the dynamical system for with jump.
Table 10 shows the estimates of the suboptimal size and the parameters ,, , , , , , , , , ,, , , and , for each of the energy commodity data sets. According to (58), the estimate , , is positive if commodity is cooperating with commodity , and negative if commodity is competing with commodity . There is no interaction between the two commodities if . It is apparent from the graph (from in Column 6) that coal is competing with natural gas during this period because the estimates of are mostly negative. It is apparent that natural gas and crude oil are either cooperating or competing, depending on the time period.
The graph of the drift coefficient’s parameters with estimates in Table 10 are given in Figure 11 below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Figures 11(a)–11(i) show the graph of the suboptimal interaction coefficient parameters , , , , , , , , and . According to (58), the estimate , , is positive if commodity is cooperating with commodity , and negative if commodity is competing with commodity . There is no interaction between the two commodities if . It is apparent from the graph of that coal is competing with natural gas because the estimates of are mostly negative. This shows that for fixed price of natural gas, for every increase in the price of coal, the price of natural gas decreases on the average by . Also, it is apparent that natural gas and crude oil are either cooperating or competing, depending on the time period.
Figures 12(a)–12(c) are the graphs of , , and against time for the daily Henry Hub natural gas price data set [36], daily crude oil price data set [37], and daily coal price data set, respectively.
(a)
(b)
(c)
Table 11 gives the suboptimal estimates of the parameters , , , , , , , , , , , and for each of the energy commodity data sets.
Figures 13(a)–13(c) are the graphs of , , , , , , , , and against time for the daily Henry Hub natural gas price data set [36], daily crude oil price data set [37], and daily coal price data set, respectively.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
Table 12 shows the real and simulated estimates for the spot price processes , .
Figure 14 shows the graph of the real and simulated prices (with jump) for natural gas, crude oil, and coal data set.
(a)
(b)
(c)
Figures 14(a)–14(c) show the graph of the real and simulated spot prices for the daily Henry Hub natural gas data set [36], daily crude oil data set [37], and daily coal data set [35], respectively. The red line represents the real data set , while the blue line represent the simulated data set . It is clear that the graph fits well. To reduce magnitude of error, we increase the magnitude of time delay. It is obvious that these curves fit better than the curves in Figure 7. It follows that the interconnected dynamical system with jump process incorporated into it performs better than the one without jump.
8. Forecasting
In this section, we shall sketch an outline about forecasting problem for the case where there is no jump. The sketch for the case where jump exist is similar. An suboptimal simulated value (, ) at time , , is used to define a forecast (, ) for () at the time for the system of energy commodity model.
8.1. Forecasting for Energy Commodity Model
In the context of Section 6.1, for , we begin forecasting from time . Using the data set up to time , we compute , , , , , , , , and for . We assume that we have no information about the real data set . Under these considerations, imitating the computational procedure outlined in Section 7 and using solutions to (66) and (67), we find the estimate of the forecast and at time as follows:
To find the forecasted value and , we need to re-estimate and by first simulating as follows:
Imitating argument in (84), we define
Thus, and . where the estimates , , , , , and , are estimated with respect to the known past data set up to the time . We note that is the suboptimal estimate for at time .
To determine (), we need , , , , , and , . Since we only have information of real data up to time , we use the forecasted estimate as the estimate of at time and to estimate , , , , , and , . Hence, we can write as
For , to find (, ), we use the estimates
Continuing this process in this manner, to find , we use the estimates
8.1.1. Prediction/Confidence Interval for Energy Commodities
In order to be able to assess the future certainty, we also discuss about the prediction/confidence interval. We define the confidence interval for the forecast of the state () at time , , as where is the estimate for the sample standard deviation for the forecasted state , and is the estimate for the sample standard deviation for the forecasted state derived from the following iterative process
It is clear that the 95% confidence interval for the forecast at time is where the lower ends denote the lower bounds of the state estimate and the upper ends denote the upper bounds of the state estimate.
8.1.2. Illustration: Prediction/Confidence Interval for Energy Commodities with No Jump
For the case of no jump, the following graphs show the simulated, forecasted, and 95 percent confidence limit for the daily Henry Hub natural gas data set [36], daily crude oil data set [37], and daily coal data set [35], respectively.
Figures 15(a)–15(c) show the graph of the forecast and 95 percent confidence limit for the case where there is no jump for the daily Henry Hub natural gas data set [36], daily crude oil data set [37], and daily coal data set [35], respectively. The figure shows two regions: the simulation region and the forecast region . For the simulation region , we plot the real data set together with the simulated data set as described in Figure 14. For forecast region , we plot the estimate of the forecast as explained in Section 8.1.1.
(a)
(b)
(c)
8.1.3. Illustration: Prediction/Confidence Interval for Energy Commodities with Jump
For the jump process, the following graphs show the forecast and 95 percent confidence limit for the daily Henry Hub natural gas data set [36], daily crude oil data set [37], and daily coal data set [35], respectively.
Figures 16(a)–16(c) show the graphs of the forecast and 95 percent confidence limit for the case where there is jump for the daily Henry Hub Natural gas data set [36], daily Crude Oil data set [37], and daily Coal data set [35], respectively. Figures 16(a)–16(c) show two regions: the simulation region and the forecast region . For the simulation region , we plot the real data set together with the simulated data set as described in Figure 14. For the forecast region , we plot the estimate of the forecast as explained in Section 8.1.1.
(a)
(b)
(c)
9. Conclusion and Future Work
We examine the relationship between the prices of the energy commodities: natural gas, crude oil, and coal. To do this, multivariate stochastic models with and without external random interventions are developed to describe the relationship between the price of energy commodities. The time-varying parameters in the stochastic model are estimated using the LLGMM method. Figures 7 and 14 suggest that the model with jump fits better than the model without jump. For , the estimates for the drift interaction parameters for the case where jump is incorporated suggest that there are definitely interactions between the spot price of these three commodities. As discussed in (57) and (58), the sign of these parameters suggest that there is competition or cooperation between commodities and . The estimate of the parameter in Tables 3 and 10 and Figures 4 and 11 suggests that these commodities either compete or cooperate with each other depending on the time period. We can also describe the relationship between any two commodity and , , based on the overall average . For example, for the case where jump is incorporated, and . This suggests that on the average, there is competition between these two commodities. Also, and . This indicates that on the average, there is cooperation between natural gas and crude oil. Finally, and . Therefore, on the average, there is competition between crude oil and coal. In the future, we plan to apply the local lagged adapted generalized method of moments to interconnected nonlinear stochastic dynamic model for log-spot price, expected log-spot price, and volatility process. Also, we plan to incorporate delay in the multivariate interconnected nonlinear stochastic model. We plan to be able to apply the extended local lagged adapted generalized method of moments to other multivariate interconnected nonlinear dynamic model different from energy commodity model.
Data Availability
We will provide the data upon request.
Disclosure
We would also like to mention that this research paper was presented in the Ph.D. Dissertation [15] and at the American Mathematical Society (AMS) conference at San Antonio, Texas.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This research is supported by the Mathematical Sciences Division, The U.S. Army Research Office Grant Number W911NF-12-1-0090.