Subordinated affine structure models for commodity future prices

Kateregga, M.; Mataramvura, S.; Taylor, D. Article Subordinated affine structure models for commodity future prices Cogent Economics & Finance Provided in Cooperation with: Taylor & Francis Group Suggested Citation: Kateregga, M.; Mataramvura, S.; Taylor, D. (2018) : Subordinated affine structure models for commodity future prices, Cogent Economics & Finance, ISSN 2332-2039, Taylor & Francis, Abingdon, Vol. 6, Iss. 1, pp. 1-26, https://doi.org/10.1080/23322039.2018.1512360 This Version is available at: https://hdl.handle.net/10419/245156 Standard-Nutzungsbedingungen: Die Dokumente auf EconStor dürfen zu eigenen wissenschaftlichen Zwecken und zum Privatgebrauch gespeichert und kopiert werden. Sie dürfen die Dokumente nicht für öffentliche oder kommerzielle Zwecke vervielfältigen, öffentlich ausstellen, öffentlich zugänglich machen, vertreiben oder anderweitig nutzen. Sofern die Verfasser die Dokumente unter Open-Content-Lizenzen (insbesondere CC-Lizenzen) zur Verfügung gestellt haben sollten, gelten abweichend von diesen Nutzungsbedingungen die in der dort genannten Lizenz gewährten Nutzungsrechte. Terms of use: Documents in EconStor may be saved and copied for your personal and scholarly purposes. You are not to copy documents for public or commercial purposes, to exhibit the documents publicly, to make them publicly available on the internet, or to distribute or otherwise use the documents in public. If the documents have been made available under an Open Content Licence (especially Creative Commons Licences), you may exercise further usage rights as specified in the indicated licence. https://creativecommons.org/licenses/by/4.0/ Full Terms & Conditions of access and use can be found at https://www.tandfonline.com/action/journalInformation?journalCode=oaef20 Cogent Economics & Finance ISSN: (Print) 2332-2039 (Online) Journal homepage: https://www.tandfonline.com/loi/oaef20 Subordinated affine structure models for commodity future prices M. Kateregga, S. Mataramvura & D. Taylor | To cite this article: M. Kateregga, S. Mataramvura & D. Taylor | (2018) Subordinated affine structure models for commodity future prices, Cogent Economics & Finance, 6:1, 1512360, DOI: 10.1080/23322039.2018.1512360 To link to this article: https://doi.org/10.1080/23322039.2018.1512360 © 2018 The Author(s). This open access article is distributed under a Creative Commons Attribution (CC-BY) 4.0 license. Published online: 07 Sep 2018. Submit your article to this journal Article views: 380 View related articles View Crossmark data FINANCIAL ECONOMICS | RESEARCH ARTICLE Subordinated affine structure models for commodity future prices M. Kateregga, S. Mataramvura and D. Taylor Cogent Economics & Finance (2018), 6: 1512360 FINANCIAL ECONOMICS | RESEARCH ARTICLE Subordinated affine structure models for commodity future prices M. Kateregga 1 *, S. Mataramvura 1 and D. Taylor 1 Abstract: To date the existence of jumps in different sectors of the financial market is certain and the commodity market is no exception. While there are various models in literature on how to capture these jumps, we restrict ourselves to using subordinated Brownian motion by an α-stable process, α∈(0,1), as the source of randomness in the spot price model to determine commodity future prices, a concept which is not new either. However, the key feature in our pricing approach is the new simple technique derived from our novel theory for subordinated affine structure models. Different from existing filtering methods for models with latent variables, we show that the commodity future price under a one factor model with a subordinated random source driver, can be expressed in terms of the subordinator which can then be reduced to the latent regression models commonly used in population dynamics with their parameters easily estimated using the expectation maximisation method. In our case, the underlying joint probability distribution is a combination of the Gaussian and stable densities. Subjects: Applied Mathematics; Financial Mathematics; Mathematical Finance; Quantitative Finance; Statistics & Probability; Statistics; Statistics for Business, Finance & Economics Keywords: stable distributions; affine structure models; latent regression models AMS Subject Classification: 62G05; 62G07; 62G32 ABOUT THE AUTHOR M. Kateregga holds a Ph.D. in Mathematical Finance from the University of Cape Town. Michael is a Software Developer at Mira Networks in South Africa. His main interests include financial markets research and functional programming particularly Haskell. PUBLIC INTEREST STATEMENT This paper is entitled Subordinated Affine Structure Models for Commodity Futures Prices. The aim is to provide a mathematical background on a special family of processes that could capture unusual patterns such as extreme events in the commodities market. Such events could range from wild fires on farms to Tsunamis and as a result the commodities spot prices would be affected significantly. This realisation provides the motivation of the study in this paper. As a risk mitigation process, it is undeniable that investors in the commodities market would be interested in how this kind of risk could be captured in the future pricing model and more so, how the model would be calibrated to historical data. The current paper seeks to answer these questions and establishing a novel approach on pricing commodities futures. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 2 of 26 Received: 02 April 2018 Accepted: 11 August 2018 First Published: 20 August 2018 *Corresponding author: M. Kateregga, University of Cape Town, Rondebosch, Cape Town, 7700, South Africa E-mail: [email protected] Reviewing editor: Lanouar Charfeddine, Qatar University, Qatar Additional information is available at the end of the article © 2018 The Author(s). This open access article is distributed under a Creative Commons Attribution (CC-BY) 4.0 license. 1. Literature A large volume of literature on commodities market has been published since the invention of the continuous benchmark model of Black and Scholes (1973) for pricing options and corporate liabilities. Among the many models developed, the widely used and referenced study on commodities is the work by Schwartz (1997). From the latter, numerous models have been developed as a result of the growing commodities market in terms of volumes traded and complexities of their contracts over the years. We give an account of the various literatures relevant to this research in Table 1. The current research builds on findings from Kyriakou, Nomikos, Pouliasis, and Papapostolou (2016) by extending the results to subordinated Brownian motion. A selection of key commodity jump models that have developed overtime. 2. Introduction Commodities exhibit distinctive features that a good model ought to capture. For instance to estimate the commodities market as closely as possible, one has to factor in jumps in the underlying spot price. However, models designed with a jump component are non-trivial. In this research, we derive a relatively easy estimation method for commodities prices using subordination as a proxy for introducing jumps. Other features include mean-reversion, contango, backwardation and seasonality. Commodities also experience extreme volatility and price spikes resulting in heavy-tailed distribution of the returns. Commodity markets are unique compared with other markets such as equity, bond, currency or interest rate markets in the sense that most commodities are real physical assets that are produced, transported, stored and consumed. They are not assets valued on long-lived companies like in equity markets. As indicated in Fama and French (1988), commodity pricing can be approached from two perspectives, the theory of storage which explains why high supplies and inventories running at minimum would result into contango, low futures and spot price volatilities, and in turn futures premiums being equivalent to full storage costs. On the other hand, why low supplies and enhanced production inventory levels yield to backwardation, a rise in volatilities of the spot and the nearby future prices. Another feature explained by this theory is the periodically continuously compounded convenience yield (usually denoted by δ) on inventory which is the benefit of holding a physical commodity as opposed to having a future contract of its delivery at some future time and second, the cost of storage. The future price motivated by the theory of storage is given by Ft;T¼SteðrcδÞðTtÞ;(1) Table 1. Various commodity jump models Author(s) (Year) Study Result G´omez-Valle and MartnezRodr`ıguez (2017) Jump size distribution More accurate than a diffusion model Kateregga, Mataramvura, and Taylor (2017) Jump distribution Better estimation method than Maximum Likelihood Kyriakou et al. (2016) Jumps & Heston-type stochastic vol. Evidence of jumps and stochastic vol. in oil Prokopczuk, Symeonidis, and Simen (2015) Characterisation of jumps Jump importance to vol. Schmitz, Wang, and Kimn (2014) Stoch. vol., price jumps, seasonality & stoch. cost of carry A model that capture all the four features. Maslyuka, Rotarub, and Dokumentovc (2013) Frequency of price discontinuities Evidence of discontinuities Deng (2000) Multiple jumps, regime–switching & stoch. vol. Prices of various commodities. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 3 of 26 where caccounts for the storage costs, ris the periodically continuously compounded interest rate, Stis the current spot price and Tis the maturity date of the future contract. The second perspective is the theory of expected risk premium discussed in Keynes (1930) and Hicks (1939). It asserts that the future prices are given by the discounted (by the risk premium) expected future spot price: Ft;T¼Et½STerγ½Tt;(2) where γis the risk premium and Et½ ¼ E½jFt,Ftis the filtration up to time t. A number of models based on the latter approach have been developed over the years to mimic the market as closely as possible for various commodities. This includes Schwartz’s common continuous stochastic factor models Schwartz (1997), Schwartz and Smith (2000) and the jump models of Kyriakou et al. (2016). The motivation and contribution of this paper are based on the existing erratic features in electricity and energy markets, where jumps are evident, resulting in skewed distributions of the spot prices. We consider a subordinated Brownian motion by an α-stable process, α{ð0;1Þ, as the source of randomness in the underlying asset to model commodity future prices. The stunning feature in our pricing approach is the new simple technique derived from our novel approach for subordinated affine structure models. We show that the affine property is attainable and applicable to generalised commodity spot models, and as illustration, we consider a stochastic differential equation with subordinated Brownian motion as the source of randomness to derive the commodity future prices. It is argued in some existing literature that the likelihood function exists in integrated form for models with singular noise meanwhile for cases of partially observed processes, a filtering technique is required, see for instance Date and Ponomareva (2010) and Yang, Lin, Luk, and Nahar (2014). However, the work presented in this paper provides a new approach of pricing commodity futures for models with latent variables using the maximum expectation maximisation. We show that the commodity future price under a one-factor model with a subordinated random source driver can be expressed in terms of the subordinator, which can then be reduced to the latent regression models commonly used in population dynamics with their parameters estimated easily using the expectation maximisation method. In our case, the underlying joint probability distribution is a combination of Gaussian and stable densities. The rest of the paper is organised as follows. The following Section 3introduces some features of stable processes essential to this work. In Section 4, we review the concept of affine models and extend the idea of obtaining Laplace transforms of random processes to subordinated processes. In Section 5, we derive our pricing formulas for commodity futures using the results derived in Section 4. In Section 6, we discuss the numerical implementation of our one-factor commodity future models. Section 7concludes. 3. Stable processes The discussion in this section is mainly based on Kateregga et al. (2017). A stable or α-stable process, α{ð0;2, belongs to the general class of Lévy distributions. It has a limiting distribution with a definitive exponent parameter αthat determines its shape. The following two definitions follow from Rachev (2003). Definition 3.1 Let X1;X2;;Xnbe independent and identically distributed random variables and suppose a random variable S defined by Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 4 of 26 1 an ∑ n i¼1 Xibn ! !S;(3) where “→”represents weak convergence in distribution, anis a positive constant and bnis real. Then, Sis a stable process. The constants anand bnneed not to be finite. Definition 3.1 allows the modelling of a number of natural phenomena beyond normality using stable distributions. The fact that anand bndo not necessarily need to be finite provides the generalised central limit theorem. Definition 3.2 (Generalised Central Limit Theorem) Suppose X1;X2, denotes a sequence of identically distributed independent random variables from some arbitrary distribution and let sequences an{Rand bn{Rþ. Then, we define a sequence Zn:¼1 bn ∑ n i¼1 Xian ! (4) of sums Znsuch that their distribution functions weakly converge to some limiting distribution. That is PðZn<xÞ!HðxÞ;n!1;(5) where H(x) is some limiting distribution. The traditional central limit theorem assumes finite mean a:¼E½Xiand finite variance σ2:¼ Var½Xiand defines the sequence of sums Zn:¼1 σffiffiffi n p∑ n i¼1 Xina ! ;(6) such that the distribution functions of Znweakly converge to hsGðxÞ. That is Pðx1<Zn<x2Þ!ðx2 x1 hsGðxÞdx;n!1 (7) where hsGðxÞdenotes the standard Gaussian density hsGðxÞ¼ 1 ffiffiffiffiffiffi 2π pexpðx2=2Þ:(8) Suppose the identically distributed independent random variables Xiequal to a positive constant c almost surely and the sequences anand bnin (4) are defined by an¼ðn1Þcand bn¼1, then Zn is also equal to cfor all n>0 almost surely. In this case, the random variables Xiare mutually independent and as a consequence, the limiting distribution for the sums Znbelongs to the stable family of distributions by definition. Therefore, stable distributions behave similarly to the central limit theorem for distributions with a finite second-order moment (the Gaussian), Crosby (2008). This is one reason why they are regarded as stable. They are also preferred compared with all other laws such as the Normal Inverse Gaussian (NIG), Variance Gamma (VG) and other distributions from the generalised hyperbolic family because of their heavier tails. Definition 3.3 Samoradnitsky and Taqqu (1994)Anα-stable distribution is a four-parameter family of distributions denoted by Sðα;β;ν;μÞwhere (1) α{ð0;2is the characteristic exponent responsible for the tail of the distribution. (2) β{½1;1is responsible for skewness. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 5 of 26 (3) ν>0 is the scale parameter (sometimes referred to as variance when α¼2). (4) μ{Ris the location (sometimes referred to as mean). Densities of α-stable distributions do not have closed-form representations except the Gaussian (α¼2), Cauchy (α¼1) and Inverse Gaussian or Pearson (α¼0:5) distributions. The analysis of stable processes is usually through their characteristic functions and Laplace or Fourier transformation. Unlike their densities, their characteristic functions always exist. Literature on their integral representations and density functions is provided in Zolotarev (1964), 1980, Zolotarev (1986)). The distribution functions for the different αvalues have been tabulated in Dumouchel (1971), Fama and Roll (1968) and Holt and Crow (1973). Definition 3.4 (Gajda and Wyłoman`Ska (2012)) Let Stand Ltdenote an α-stable process and its respective inverse. Then for t{½0;T, we define the process Ltas 1 Ls:¼inf t:St>s fg if s{½0;StÞ Tif s¼ST: (9) Stand Ltare non-decreasing and cádlág with their graphical representations given in Figure 1. 3.1. Density and characteristic functions Let ðXt;t0Þdenote a Lévy process. The characterisation of Xtis usually given by the Lévy– Khintchine formula. Definition 3.5 (Lévy–Khintchine) Applebaum (2004) Consider a Lévy process X¼ðXtÞt0. There exists b{Rand σ0 such that the characteristic function of Xis ΦðtÞ:¼E½eitX¼exp itb 1 2σ2t2þðR0 fg ðeitx 1itxIxjj<1ÞðdxÞ  ;(10) where Iis the indicator function and is a σ-finite measure satisfying the constraint t 0 0.2 0.4 0.6 0.8 1 1.2 α-Stable Process 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 t 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.2 0.4 0.6 0.8 LtSt 1 1.2 1.4 Inverse α-Stable Process Figure 1. The top graph shows a plot of a stable process Stand the bottom graph shows its inverse process Ltsimulated using exponent parameter value, α¼0:8, plotted against time on the horizontal. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 6 of 26 ðRd0 fg minð1;jxj2ÞðdxÞ<1;alternatively ðRd0 fg jxj2 1þjxj2ðdxÞ<1:(11) Definition 3.6 (The Lévy-Itô Decomposition) Applebaum (2004)IfXtis a Lévy process, there exists b{R, a Brownian motion BσðtÞwith variance σ{Rþand an independent Poisson random measure N on RþðR0 fgÞsuch that, for each t0, Xt¼bt þBσðtÞþðjxj<1 x~ Nðt;dxÞþðjxj1 xNðt;dxÞ;(12) where b¼EX1ðx jj 1 xNð1;dxÞ  :(13) To preserve the martingale property, the compound Poisson random measure is compensated as ~ N¼Ntλwhere λis a Lévy measure satisfying (11). For process St, we require σ¼0 in (10) or B¼0 in (12) and the Lévy measure in (11) given by ðxÞ¼ C jxj1þαdx;C>0;(14) The characteristic function Φof Stis obtained using the domain of attraction of stable random variables (See Grigelionis, Kubilius, Paulauskas, Statulevicius, & Pragarauskas, 1999) and the Lévy– Khinchine representation formula (See Definition 3.5 or Applebaum (2004) for a detailed explanation), i.e. ΦðθÞ¼E½expðiθSÞ ¼ expðναjθjα½1iβsignðθÞtanðπα 2ÞþiμθÞ;for αÞ1: expðνθ jj½ 1þiβsignðθÞ2 πlog θ jjþiμθÞ;for α¼1: (15) Using Fourier transformation, the density function of Stis given by hSðt;uÞ¼1 πð1 0 eiuθΦðθÞdθ:(16) From Definition 3.4, it is easy to see the equivalence relation Su<t,Ltu:(17) It follows that Fðt;uÞ:¼PðSu<tÞ¼PðLtuÞ¼ð1 0 hLðτ;tÞτwhere hLðu;tÞdenotes the probability density function of Lt. Consequently hLðu;tÞ¼@Fðt;uÞ @u¼@ @uðt 1 hSðτ;uÞdτ:(18) According to Meerschaert and Straka (2013), the density hðt;uÞcan also be given by hSðt;uÞ¼u1=αhðtu1=αÞ;(19) where hðτÞis the density of a standard stable process with a Laplace transform ~ hðτÞ¼expðταÞ. This follows from the fact that Suhas the same distribution as u1=αS1.Asa result, the density of the inverse stable process Ltcan be given in terms of the standard stable process by hLðu;tÞ¼t αu11=αhðtu1=αÞ:(20) Figure 2shows density graphs of Stfor different exponent parameter values. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 7 of 26 6. Numerical implementation We focus on the one-factor model to explain our approach for estimating the model parameters. The data used in this section are obtained from the US Energy Information Administration and include future prices of Crude Oil (Light-Sweet, Cushing, Oklahoma) from 30 March 1983 to 6 December 2016 (8452 observations), Reformulated Regular Gasoline (New York Harbor) from 3 December 1984 to 31 October 2006 (5492 observations), Heating Oil (New York Harbor) from 2 January 1980 to 6 December 2016 (9262 observations) and Propane (Mont Belvieu, Texas) from 17 December 1993 to 18 September 2009 (3941 observations). The parameters in the seasonality function (56) are estimated by fitting the function to the historical spot prices. The spot prices used include Crude Oil from 2 January 1986 to 12 December 2016 (7807 observations), RBOB Regular Gasoline from 11 March 2003 to 12 December 2016 (3460 observations), No. 2 Heating Oil from 2 June 1986 to 12 December 2016 (7683 observations) and Propane from 9 July 1992 to 12 December 2016 (6133 observations). fðtÞ¼δ0tþδ1sinðδ2½tþ2π=264Þþδ3sinðδ4½tþ4π=264Þ: The accuracy of the fitting in Figure 3depends on the choice of the initial parameters δ0,δ1,δ2,δ3 and δ4. Figure 3. Seasonality is captured by the function fðtÞ defined in the following table. The best fit of fðtÞcan be obtained by obtaining an optimal set of the δparameters. Table 2. Estimation of parameters in the seasonality function Commodity δ0δ1δ2δ3δ4 Crude Oil 0.01097129 −0.04210297 13.29908652 18.43819425 6.28446063 Gasoline 0.001018775 −0.515953818 6.27995144 1.252457783 6.28435941 Heating Oil 0.000287309 0.159154792 12.57029401 −0.49422029 6.28182486 Propane 0.000243901 −0.252880065 12.56438012 0.071328671 6.286055367 Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 14 of 26 6.1. Equivalent latent regression model The de-seasonalised one-factor future price given by (58) can be written as y¼aþbx þε;(63) where y¼ln F,x¼Stand aand bare given by a¼ 1 2σ2þβ  α þβα  γð1eκTÞþ1 4κσ21 2σ2þβ  α þβα  2 ð1e2κTÞ(64) b¼ 1 2σ2þβ  α þβα  eκT;(65) and εis an independent random error distributed as Nð0;ΘÞwith zero mean and variance Θ. Clearly, (63) belongs to the class of latent regression models since x¼Stis not observable. There is literature where this kind of problem is handled using expectation-maximisation (EM) algorithms (see Dempster, Laird, and Rubin (1977)) to estimate the model parameters. The latent variable xcan be considered binary where in this case the EM algorithm would give estimates for a two-component normal mixture model. On the other hand, xcan be allowed to be continuous between 0 and 1 with a beta distribution as in Tarpey and Petkova (2010). The EM algorithm for estimating the model parameters in this case is more involved than for the two-component mixture model and more computationally challenging, but can be done nonetheless. Basically, the latent or unobserved xvariables are imputed by their conditional expectation given the outcomes y.Our approach is adaptable to the latter approach through the Dynkin–Lamperti Theorem (see Gupta and Nadarajah (2004)) where the unobserved variable follows a stable distribution defined on ð1;1Þ with α{ð0;2and the observable variable yrepresents the log-returns of the future prices. The Figure 4. Detection of jumps in crude oil future prices. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 15 of 26 algorithm is applied to the joint likelihood of the response y. We assume the error εis independent of the latent predictor x. The joint density for xand yis given by hðx;yÞ¼hðyjx;a;b;ΘÞSðx;α;β;ν;μÞ: ¼Nðy;aþbx;ΘÞSðx;α;β;ν;μÞ;(66) where Sðx;α;β;ν;μÞis the αstable distribution, Nðy;A;BÞdenotes the normal distribution of a random variable ywith mean Aand variance B; and Θis the variance of the outcome sample data. The marginal density of the response yis fðyÞ¼ð1 ffiffiffiffiffiffi 2π pΘ1=2expðΘ1ðyβ0β1xÞ0Θ1ðyβ0β1xÞ=2Þ Sðx;α;β;ν;μÞdx:(67) Density (67) is an example of infinite mixture models used in ecological statistics (Fox, NegreteYankelevich, and Sosa (2015)). 6.2. The EM algorithm For a data set ðx1;y1Þ;;ðxn;ynÞin (63). The log-likelihood is derived from (66) as Lðα;β;ν;μ;Θ;a;bÞ¼n 2logð2πÞn 2log Θ jj ∑ n i¼1ðyiβ0β1xiÞ0ðyiβ0β1xiÞ=2 þ∑ n i¼1 logðSðxi;α;β;ν;μÞÞ: (68) Since xis not observable, the EM algorithm requires maximising the conditional expectation of the log-likelihood given the response vector y. That is E½Lðα;β;ν;μ;Θ;a;bÞjy;(69) where at each iteration of the EM algorithm, the above conditional expectation is computed using the current parameter estimates. This current expectation-maximisation problem is similar to the problem handled in Tarpey and Petkova (2010), which implies that a similar EM algorithm can be applied here. That is, suppose ρ¼a b  is a 2 pdimension coefficient matrix where each of the p columns of ρprovides the intercept and slope regression coefficients for each of the ρresponse variables. Then the design matrix is denoted by Xwhose first column consists of ones for the intercept and the second column consists of the latent predictors xi;i¼1;;n. The multivariate regression model follows: Y¼Xρþε:(70) Moreover, and as indicated in Tarpey and Petkova (2010), the likelihood for the multivariate normal regression model can be given as Lðρ;ΘÞ¼ð2πÞnp=2jΘj1=2expðtr½Θ1ðYXρÞ0ðYXρÞ=2Þ:(71) The EM approach requires that we maximise the expectation of the logarithm of (71) conditional on Ywith respect to ρand Θ. This leads to the following optimal factors: ^ ρ¼ð~ X0XÞ1~ X0Y:(72) ^ Θ¼Y0Y^ ρ0ð~ X0XÞ^ ρ;(73) where ~ X¼E½XjYand ð~ X0XÞ¼E½X0XjY. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 16 of 26 To implement the EM method in the R programming language, we first highlight that there are a minor differences to bear in mind before implementing the algorithm as we explain in the following. First and foremost, the density of the predicator in our case is from a stable distribution. Recall, in general, the densities of stable processes cannot be expressed analytically, which makes it difficult to compute the log-likelihood. However, with the help of inbuilt packages in R including stabledist and StableEstim, the log likelihood can be satisfactorily estimated using EstimðÞ to obtain the stable parameters of St, and dstableðÞ for its corresponding stable density. Second, from the log-likelihood expression, we notice that we only require estimates of the conditional expectations of x,x2and log xwith respect to the joint probability density given the response vector y. On the other hand, we retain some of the steps in Tarpey and Petkova (2010). The initial values for the regression parameters aand bcan be obtained from fitting a two-component finite mixture model or by a preliminary search over the parameter space. Initial values for Θcan be obtained using the sample covariance matrix from the raw data. 6.3. Data for the EM algorithm The data used is stored in a data frame with three columns containing futures log-returns, spot price log-returns and binary data of 1’s and 0’s representing whether or not a jump has occurred within a given window size (see Table 4). Table 2shows the estimated parameters as a result. The parameters were obtained from 5000 data points of crude oil log-returns arranged as in Table 4. We have displayed results from only two iterations because for large data sets the code tends to be slow in addition to suffering convergence issues. However, this can be improved and using faster machines. The jump occurrence due to Stis determined by the method discussed in Lee and Mykland (2008) (also see Maslyuka et al., 2013). That is the realised return at any given time is compared with a continuously estimated instantaneous volatility σtito measure local variation arising from the continuous part of the process. The volatility σtiis estimated using a modified version of realised bipower variation calculated as the sum of products of consecutive absolute returns in the local window (see Barndorff-Nielsen & Shephard, 2004). Then, the jump detection statistic Li{ testing for jumps in returns occurring at a time tiwithin a window size Kis calculated as the ratio of realised returns to estimated instantaneous volatility: Table 3. Parameters obtained from maximum likelihood method Parameter Estimation Parameter 1st Iteration 2nd Iteration loglike 30687 31133:2 a7:74707e05 1:17456e05 b0:00218499 0:00166723 α1:6605 1:6605 β0:0651915 0:0651915 ν0:00576286 0:00576286 μ0:000415904 0:000415904 Θ0:0106344 0:0110766 T5000 days 5000 days Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 17 of 26 Li;log Yti=Yti1 ^σti ;(74) where Ytat t0 represents the commodity spot price and ^σtiis estimated by ^σ2 ti ;1 K2∑ i1 j¼iKþ2 logðYtj=Ytj1Þ logðYtj1=Ytj2Þ :(75) Care must be taken in choosing K, it must be large enough to accurately estimate integrated volatility but small enough for the variance to be approximately constant. In other words, K should be large enough but smaller than N, the number of observations so that the effect of jumps on estimating instantaneous volatility disappears. Some authors recommend Kto be computed as K¼ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 252 n p, where nis the daily number of observations, whereas 252 is the number of days in the (financial) year. Moreover, the window size should be such that K¼ OðΔtλÞwith 1<λ<0:5. For high-frequency data, Lee and Mykland (2008) recommend, for returns sampled at frequencies of 60, 30, 15 and 5 minutes, the corresponding values of Kto be 78, 110, 156 and 270. For our case, we shall choose K¼4 for crude oil future prices with returns sampled daily. 6.3.1. Detection of jumps in the data The test statistic Lfollows approximately a normal distribution when the data set has no jumps and its value becomes large otherwise. According to Lee and Mykland (2008), the region for Lis chosen based on the distribution of its maximum. For instance, suppose a particular interval ðti1;tihas no jumps and the distance between two consecutive observations in this interval is small (i.e. Δ!0), then the maximum should converge to the Gumbel variable: maxi{ ANLi jjcN sN!; (76) where has a cumulative distribution function PðxÞ¼expðexÞ, ANis the set of i{1;2;;N fg such that there is no jump in ðti1;tiand cN,sNare defined as cN¼ð2 log NÞ1=2 0:8log πþlogðlog NÞ 1:6ð2 log NÞ1=2:(77) Table 4. A snapshot of the structure of the data used Futures log-returns Spot log-returns Jump detection −0.015646306 −0.020072054 0 0.0322302 0.047966668 0 0−0.016401471 0 −0.001778683 0.001489862 0 −0.017271048 −0.019778648 1 −0.02611648 −0.022030684 1 −0.043103026 −0.033698961 0 0−0.010387745 0 0.001808809 0.003969483 0 0.037999099 0.036829588 1 −0.030483308 −0.030057522 0 0.024152859 0.019704761 0 0.009625568 0.008368168 0 Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 18 of 26 sN¼1 0:8ð2 log NÞ1=2:(78) The test is conducted by comparing the standardised maximum of Liin (76) to the critical values from the Gumbel distribution where the null hypothesis of no jump is rejected when the jump statistic is given by Li>G1ð1λÞsNþcN;(79) where G1ð1λÞis the ð1λÞquantile function of the standard Gumbel distribution. Suppose λ¼0:1, then we reject the null hypothesis of no jump when Li>sNηþcNwhere ηis such that expðeηÞ¼1η¼0:9. That is η¼logðlogð0:9ÞÞ¼2:25. Figure 5shows a graph of jumps detected in crude oil future prices where we have used 1’s to record a jump occurrence and 0’s for no jump. 7. Summary We have shown that the affine property is attainable and applicable to generalised spot models. We considered a stochastic differential equation with the source of randomness as subordinated Brownian motion as a specific example to derive the futures price. Moreover, it has been argued in some existing literature that the likelihood function exists in integrated form for models with singular noise meanwhile for cases of partially observed processes, a filtering technique is required. However, the work presented in this paper provided a new approach of pricing commodity futures for models with latent variables using the maximum expectation-maximisation, without using any filtering. Our approach is easy to implement once the joint probability density is established. The numerical implementation of the two factor model is left for future work. Acknowledgements This work was supported by funds from the National Research Foundation of South Africa (NRF), the African Institute for Mathematical Sciences (AIMS) and the African Collaboration for Quantitative Finance and Risk Research (ACQuFRR), which is the research section of the African Institute of Financial Markets and Risk Management (AIFMRM), which delivers postgraduate education and training in financial markets, risk management and quantitative finance at the University of Cape Town in South Africa. All the numerical results were generated using MATLAB and the Robust Analysis software package in R accessible at www.RobustAnalysis.com. Funding This work was supported by the University of Cape Town [OAPF-713]; Author details M. Kateregga 1 Figure 5. The large spikes in the jump statistic graph reflect the extreme events in the price where as the blue strips represent the smaller jumps that might go unnoticed. The former jumps are easily detected in returns yet the latter are not that visible. Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 19 of 26 E-mail: [email protected] ORCID ID: http://orcid.org/0000-0003-1992-1929 S. Mataramvura 1 E-mail: [email protected] D. Taylor 1 1 Software Development, University of Cape Town, Cape Town, 7700, South Africa. Citation information Cite this article as: Subordinated affine structure models for commodity future prices, M. Kateregga, S. Mataramvura & D. Taylor, Cogent Economics & Finance (2018), 6: 1512360. Note 1. The process Lis also interpreted as the first passage time of the stable process S. Correction This article has been republished with minor changes. These changes do not impact the academic content of the article. References Applebaum, D. (2004). L´evy processes and stochastic calculus. Cambridge studies in advanced mathematics. Cambridge: Cambridge University Press. Barndorff-Nielsen, O. E., & Shephard, N. (2004). Power and bipower variation with stochastic volatility and jumps. Journal of Financial Econometrics,2(1), 1–37. doi:10.1093/jjfinec/nbh001 Black, F., & Scholes, M. (1973). The pricing of options and corporate liabilities. Journal of Political Economy,81(3), 637–654. doi:10.1086/260062 Bochner, S. (2012). Harmonic analysis and the theory of probability. Dover books on mathematics. Mineola, NY: Dover Publications. Crosby, J. (2008). A multi-factor jump-diffusion model for commodities. Journal of Quantitative Finance,8 (2), 181–200. doi:10.1080/14697680701253021 Date, P., & Ponomareva, K. (2010). Linear and nonlinear filtering in mathematical finance: A review. IMA Journal of Management Mathematics,22(3), 195–211. doi:10.1093/imaman/dpq008 Dempster, A. P., Laird, N. M., & Rubin, D. B. (1977). Maximum likelihood from incomplete data via the EM algorithm. Journal of the American Statistical Association,39,1–38. Deng, S. (2000). "Stochastic models of energy commodity prices and their applications: Mean-reversion with jumps and spikes," Working Paper PWP-073, University of California Energy Institute. Duffie, D., Filipovi´C, D., & Schachermayer, W. (2003). Affine processes and applications in finance. Annals of Applied Probability,13, 984–1053. doi:10.1214/ aoap/1060202833 Dumouchel, W. H. (1971). Stable distributions in statistical inference. The Journal of the American Statistical Association,78(342), 469–477. Fama, E. F., & French, K. R. (1988). Permanent and temporary components of stock prices. The Journal of Political Economy,96(2), 246–273. doi:10.1086/261535 Fama,E.F.,&Roll,R.(1968). Some properties of symmetric stable distributions. Journal of the American Statistical Association,63(323), 817–836. Fox, G., Negrete-Yankelevich, S., & Sosa, V. (2015). Ecological statistics: Contemporary theory and application. Oxford: Oxford University Press. Fusai, G., & Roncoroni, A. (2007). Implementing models in quantitative finance: Methods and cases. Springer finance. Berlin, Heidelberg: Springer. G´omez-Valle, L., Habibilashkary, Z., & Martnez-Rodr`ıguez, J. (2017). A new technique to estimate the risk-neutral processes jump-diffusion commodity futures models. Journal of Computational and Applied Mathematics, 309,435–441. doi:10.1016/j.cam.2015.12.028 Gajda, J., & Wyłoman`Ska, A. (2012). Geometric Brownian motion with tempered stable waiting times. Journal of Statistical Physics 1,148(2), 296–305. doi:10.1007/ s10955-012-0537-3 Grigelionis, B., Kubilius, J., Paulauskas, V., Statulevicius, V., & Pragarauskas, H. (1999). Probability theory and mathematical statistics: Proceedings of the seventh vilnius conference (1998), Vilnius, Lithuania, 12–18 August, 1998. Probability theory and mathematical statistics. TEV. Gupta, A., & Nadarajah, S. (2004). Handbook of beta distribution and its applications. Statistics: A series of textbooks and monographs. Abingdon: Taylor & Francis. Hicks, J. R. (1939). Value and Capital. Cambridge: Oxford University Press. Holt,D.R.,&Crow,E.L.(1973). Tables and graphs of the stable probability density function. Journal of Research of the National Bureau of Standards, Section B,77, 143–198. doi:10.6028/jres.077B.017 Kateregga, M., Mataramvura, S., & Taylor, D. (2017). Parameter estimation for stable distributions with application to commodity futures log-returns. Cogent Economics & Finance,5(1), 1318813. doi:10.1080/ 23322039.2017.1318813 Keller-Ressel, M. (2008). Affine processes theory and applications in finance. (PhD thesis), Technischen Universita¨t Wien, Vienna, Austria. Keller-Ressel, M., Schachermayer, W., & Teichmann, J. (2011). Affine processes are regular. Probability Theory and Related Fields,151(3), 591–611. doi:10.1007/s00440-010-0309-4 Keynes, J. M. (1930). A treatise on money, vol II. London, United Kingdom: Macmillan. Kyriakou, I., Nomikos, N. K., Pouliasis, P. K., & Papapostolou, N. C. (2016). Affine-structure models and the pricing of energy commodity derivatives. European Financial Management,22(5), 853–881. doi:10.1111/eufm.12071 Lee, S., & Mykland, P. A. (2008). Jumps in financial markets: A new nonparametric test and jump clustering. Review of Financial Studies,21(6), 2535–2563. Maslyuka, S., Rotarub, K., & Dokumentovc, A. (2013). Price discontinuities in energy spot and futures prices. Monash Economics Working Papers 33-13, Monash University, Department of Economics. Meerschaert, M. M., & Straka, P. (2013). Inverse stable subordinators. Mathematical Modelling of Natural Phenomena,8(2), 1–16. doi:10.1051/mmnp/20138201 Miller, J. F. (1951). Moment-generating functions and laplace transforms. Journal of the Arkansas Academy of Science,4(1), 16. Prokopczuk, M., Symeonidis, L., & Simen, C. W. (2015). Do jumps matter for volatility forecasting? Evidence from energy markets. Journal of Futures Markets,20,1–15. Rachev, S. (2003). Handbook of heavy tailed distributions in finance: Handbooks in finance. Handbooks in Finance. Amsterdam: Elsevier Science. Rouah, F., & Heston, S. (2015). The heston model and its extensions in VBA. Wiley Finance. Hoboken, NJ: Wiley. Samoradnitsky, G., & Taqqu, M. (1994). Stable non-gaussian random processes: Stochastic models with infinite variance. Stochastic modeling series. United Kingdom: Taylor & Francis. Schmitz, A., Wang, Z., & Kimn, J.-H. (2014). A jump diffusion model for agricultural commodities with Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 20 of 26 bayesian analysis. Journal of Futures Markets,34(3), 235–260. doi:10.1002/fut.21597 Schwartz,E.S.(1997). The stochastic behaviour of commodity prices: Implications for valuation and hedging. The Journal of Finance,52(3), 923–973. doi:10.1111/j.1540-6261.1997.tb02721.x Schwartz, E. S., & Smith, J. E. (2000). Short-term variations and long-term dynamics in commodity prices. Management Science,46(7), 893–911. doi:10.1287/ mnsc.46.7.893.12034 Tarpey, T., & Petkova, E. (2010). Latent regression analysis. Statistical Modelling,10(2), 133–158. doi:10.1177/ 1471082X0801000202 Yang,J.,Lin,B.,Luk,W.,&Nahar,T.(2014). Particle filtering-based maximum likelihood estimation for financial parameter estimation. In 2014 24th International Conference on Field Programmable Logic and Applications (FPL), 1–4. Zolotarev, V. (1986). One-dimensional stable distributions. Translations of mathematical mono-graphs. Providence, Rhode Island, United States: American Mathematical Society. Zolotarev, V. M. (1964). On the representation of stable laws by integrals. Trudy Matematicheskogo Instituta Imeni V.A. Steklova,71,46–50. Zolotarev, V. M. (1980). Statistical estimates of the parameters of stable laws. Banach Center Publications,6(1), 359–376. doi:10.4064/-6-1-359376 Appendix A Proof for Theorem 5.1 By applying Itô’s formula to Yt¼ln zt,it is readily seen that dYt¼κðγYtÞdStþσdBSt;(80) where γ:¼θσ2 2κand the future price with maturity date Tis given by (see (see Schwartz (1997)) FðT;zTÞ¼E½zT¼E½eYT:(81) Theorem 4.5 suggests an explicit representation of (81) is attainable and it can be deduced by considering first, the continuous case E½eXt. Suppose a continuous mean-reverting model given by dXt¼κðγXtÞdt þσdBt;X0¼x:(82) The corresponding affine forms of the coefficients according to Theorem 4.5 yield: K0¼κγ;K1¼κ;K2¼0;H0¼H1¼0;H2¼σ:(83) Since Xtis regular affine and σis constant for all t{½0;T, then we have E½eiuXt¼expðψ0ðt;uÞþψ1ðt;uÞxþψ2ðt;uÞσÞ;u{iR ¼C;(84) where the functions ψ0ðt;uÞ,ψ1ðt;uÞand ψ2ðt;uÞsatisfy the set of Riccati equations: @ψ0 @t¼KT 0ηþ1 2ηTH0η¼κγψ1;ψ0ð0;uÞ¼0:(85) @ψ1 @t¼KT 1ηþ1 2ηTH1η¼κψ1;ψ1ð0;uÞ¼iu:(86) @ψ2 @t¼KT 2ηþ1 2ηTH2η¼1 2σψ2 1;ψ2ð0;uÞ¼0:(87) where ηT¼ðψ1ψ2Þ. The solution set to the system of Riccati equations is given by ψ1ðt;uÞ¼iueκt;(88) ψ0ðt;uÞ¼iuγð1eκtÞ:(89) ψ2ðt;uÞ¼σu2 4κð1e2κtÞ:(90) Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 21 of 26 Using (84) where u¼i, one can easily deduce E½eXτleading to the price of a commodity future price under a continuous model framework. Capitalising on the affine nature of Ytand Theorem 4.6, we deduce the representation for E½eYτas: E½eiuYt¼expðψ0ðt;mðuÞÞþψ1ðt;mðuÞÞStþψ2ðt;mðuÞÞσÞ;(91) where the volatility σis a constant and the system of Riccati equations takes the form @ψ0 @tðt;mðuÞÞ ¼ κγψ1;ψ0ð0;mðuÞÞ ¼ 0:(92) @ψ1 @tðt;mðuÞÞ¼κψ1;ψ1ð0;mðuÞÞ ¼ mðiuÞ:(93) @ψ2 @tðt;mðuÞÞ ¼ 1 2σψ2 1;ψ2ð0;mðuÞÞ ¼ 0:(94) Consequently, the solution set is directly deduced from (88) (90) to obtain ψ1ðt;mðuÞÞ ¼ mðiuÞeκt;(95) ψ0ðt;mðuÞÞ ¼ mðiuÞγð1eκtÞ;(96) ψ2ðt;mðuÞÞ ¼ 1 4κσmðiuÞ2ð1e2κtÞ:(97) Setting u:¼iyields E½eYt¼expðmð1Þγð1eκtÞþ1 4κσ2mð1Þ2ð1e2κtÞþmð1ÞSteκtÞ:(98) The required result follows by substituting the Lévy exponent mð1Þfrom (41). Justifications for theorem 5.3 ψ1ðt;u1;u2Þ¼mðiu1Þeκt:(99) ϕ2ðt;u1;u2Þ¼2κmðiu1Þ υ2eκt∑ 1 j¼1 djmðiu1ÞjejκtþIfðt;u1Þ Cðu1;u2Þ1 2υ2ðt 0 Ifðs;u1Þds :(100) ϕ0ðt;u1;u2Þ¼θmðiu1Þð1eκtÞ þ2κλεmðiu1Þ υ2∑ 1 j¼1 djmðiu1Þj1 κð1þjÞð1eκtð1þjÞÞ  :(101) þðt 0 Ifðs;u1Þ Cðu1;u2Þ1 2υ2ðs 0 Ifðτ;u1Þdτ ds;(102) where coefficients dj  1 j¼1satisfy Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 22 of 26 djþ1¼ ∑ j1 i¼1 djdj11j>1ρυ κdj1j>0þυ2 4κ21j¼0 ðjþ1κλ κÞ:(103) Note that to ensure that the futures price is positive real, the values of ui{C;i¼1;2 fg have to be chosen carefully which in this case it could be ui¼i. The factor Cðu1;u2Þis defined as Cðu1;u2Þ¼ exp ρυmðiu1Þ κþ2∑ 1 j¼1 di 1þjmðiu2Þjþ1 ! mðiu2Þþ2κmðiu1Þ υ2∑ 1 j¼1 djmðiu1Þj:(104) Finally, the integrating factor Ifis such that Ifðt;u1Þ¼exp λtρυmðiu1Þ κeκtþ2eκt∑ 1 j¼1 dj jþ1mðiu1Þjþ1 ! :(105) ðt 0 Ifðτ;u1Þτ¼Ifðt;uiÞ ðλþρυmðiu1Þeκt2eκt∑ 1 j¼1 dj jþ1mðiu1Þjþ1  exp ρυmðiu1Þ κþ2∑ 1 j¼1 dj jþ1mðiu1Þjþ1 ! ðλþρυmðiu1Þ2∑ 1 j¼1 dj jþ1mðiu1Þjþ1Þ :(106) We provide the proof using the following proposition and subsequent lemmas. Proposition 8.2 The mean and variance of the model (60)–(61) are given by μ:¼κðθXtÞ λðεVtÞ  ;σ:¼ffiffiffiffiffi Vt p0 ρυ ffiffiffiffiffi Vt pυffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1ρ2 pffiffiffiffiffi Vt p  ; H:¼σσT¼VtρυVt ρυVtυ2Vt  : (107) Moreover, their affine forms can be given as linear models of both Xand V: μ¼K0þK1XtþK2Vt;(108) H¼H0þH1XtþH2Vt;(109) where K0¼κθ λε  ;K1¼κ 0  ;K2¼0 λ  :(110) H0¼00 00  ;H1¼00 00  ;H2¼1ρυ ρυ υ2  :(111) As a consequence, we deduce the following system of Riccati equations: Kateregga et al., Cogent Economics & Finance (2018), 6: 1512360 https://doi.org/10.1080/23322039.2018.1512360 Page 23 of 26