Sequential Monte Carlo Filters with Parameters Learning for Commodity Pricing Models

In this article, an estimation methodology based on the sequential Monte Carlo algorithm is proposed, that jointly estimate the states and parameters, the relationship between the prices of futures contracts and the spot prices of primary products is determined, the evolution of prices and the volatility of the historical data of the primary market (Gold and Soybean) are analyzed. Two stochastic models for an estimate the states and parameters are considered, the parameters and states describe physical measure (associated with the price) and risk-neutral measure (associated with the markets to futures), the price dynamics in the short-term through the reversion to the mean and volatility are determined, while that in the long term through markets to futures. Other characteristics such as seasonal patterns, price spikes, marketdependent volatilities, and non-seasonality can also be observed. In the methodology, a parameter learning algorithm is used, specifically, three algorithms are proposed, that is the sequential Monte Carlo estimation (SMC) for state space models with unknown parameters: the first method is considered a particle filter that is based on the sampling algorithm of sequential importance with resampling (SISR). The second implemented method is the Storvik algorithm [19], the states and parameters of the posterior distribution are estimated that have supported in low-dimensional spaces, a sufficient statistics from the sample of the filtered distribution is considered. The third method is (PLS) Carvalho’s Particle Learning and Smoothing algorithm [31]. The cash prices of the contracts with future delivery dates are analyzed. The results indicate postponement of payment, the future prices on different maturity dates with the spot price are highly correlated. Likewise, the contracts with a delivery date for the last periods of the year 2017, the spot price lower than the prices of the contracts with expiration date for 12 and 24 months is found, opposite occurs in the contracts with expiration date for 1 and 6 months.


Introduction
The models of stochastic processes for the prices of commodities with another class of market assets are differenced. The commodity market with respect to stocks, interest rates, or currencies is differenced. The financial modeling is unique, the assets are physically real, those are produced, transported, stored and consumed. In the markets of financial security are differently treated. In the big industries and the governmental organisms, the interest in the physical liquidation of the contracts, and in the delivery is focused. For this reason, a substantial amount of merchandise within and between countries must be mobilized, that a delivery cost and delays are considered, so practically impossible to liquidate these contracts in real-time, so spot contracts do not exist.

696
ESTIMATION COMMODITY PRICING MODELS the Heston model, the relationship between futures prices and spot commodity prices is proposed. An estimation methodology is developed, an estimation of states and parameters in commodity price models with structures nonlinear is proposed, that analytically in time real cannot be evaluated. An estimate in real-time by an algorithm sequential is proposed, that does not require storage of all the past data information. Interpreting various sources of uncertainty, including input, output and parameter uncertainty is realized. An approach dynamic, predict future price behavior, analyze changes in economic series structures, determines storage capacity for raw material products, are developed. The rest of the article is organized, as follows: In Section 2, the general problem is formulated; In Section 3, stochastic models of commodity prices are defined; In Section 4, a short review of particle filtering and three smoothing algorithms are proposed; In Section 5, the results are discussed; and in Section 6, a discussion and conclusions are established.

Formulation of the Problem
The financial market consists of negotiable instruments or assets, such as: stocks, bonds, currencies, commodities. Financial markets exist for the flow of money. The models to adjust the prices of primary products (agriculture, oil, copper) are different, those with respect to traditional financial models. Financial assets respect to stocks, interest rates or currencies are different, that is physically produced, transported, stored, and consumed, then, financial assets be differently treated. The state-space models in the literature are referenced, that a better performance is demonstrated, a relationship between the vector of observable time series (future prices for the delivery of the commodities) and the vector of the unobservable (called state vector) is established. Spot prices as a state variable can be viewed, which is not directly observed. Spot prices by one version discrete temporal of the stochastic process are generated, that a transition equation is called. For example, consider a general state-space model defined at discrete times, y t |x t ∼ p (y t |x t ) , observation density x t |x t−1 ∼ p (x t |x t−1 ) , state evolution density x 0 ∼ p (x 0 ) , initial state where x t are unobserved states of the system, y t are observations over time interval t ∈ {1, 2 . . . , T }, and p(x t |x t−1 ) and p(y t |x t ) state evolution and observation densities are prespecified, that can be non-Gaussian and involve nonlinearity. The state model (1) is just conditional by the previous state, and the observation density (1) is independent of previous observations and that is just conditional by the state x t at time t, y 1:t = (y 1 , . . . , y t ), and x 0:t = (x 0 , . . . , x t ). The joint distribution of states and observations by the probability chain rule are obtained, where θ known, The equation (1), that a close form has not defined, that is analytically intractable. Stochastic simulation algorithms to obtain a approximation are required. The marginal likelihood of the observations is obtained In state-space inference problems, the main interest is sequential estimation of the filtering distribution p(x t |y 1:t ).
Updating the filtering density by the standard filtering recursions can be obtained where p (y t+1 |y 1: Now, if θ and x 0:t are unknown, the evaluation of the a posteriori joint probability density p(x t , θ|y 1:t ) by Bayesian estimation techniques recursively can be obtained, the first prediction and, then updating the prediction to obtain the observation y t p(x t+1 , θ|y 1:t+1 ) = p (y t+1 |x t+1 , θ) p (x t+1 , θ|y 1:t ) p (y t+1 |y 1:t , θ) where p (y t+1 |y 1: The state-space models are appropriated, that if exists unobservable state variables, the system by a Markov process is generated. A sequential algorithm to estimate the parameters and unobservable state of the model can be used. In society, the uncertainty about volatility and the prediction of commodity prices has highly influenced, that in the agricultural, metal, energy industries, central banks, government entities, among others. In most commodity trading contracts, future prices are shown, that in cash is not available, how has to recover the true value of price, that problem is proposed. Spot prices for the pricing of derivatives and the valuation of investments are required.
To deal with the problem, a model by a system of stochastic differential equations is defined.
where • x t it is a vector of unknown states.
• θ it is a vector of unknown parameters.
is a vector of a Brownian movement.
• H k is a matrix that transforms observations into the same space of observations.
To guarantee the existence and uniqueness of the solutions of the stochastic differential equation (SDE), the fulfillment of coefficients with the global conditions and the linear growth of Lipschitz is assumed ( [16]). The particular problem, the stochastic behavior of the models are compared • The one-factor Schwartz model (1997). • The Heston model (1993).
In the first model, the logarithm of the cash prices of the merchandise by a process of reversing the Ornstein-Uhlenbeck average (O-U) is followed. The second model, the evolution of spot prices and their variance through a system of stochastic differential equations is described; also, that has nonzero price-volatility correlation, and that is essential if, the model to the skewed and leptokurtic price densities of commodity futures are captured.

Commodity price models
In this section, the commodity price models for basic products are defined, additionally, concepts for the construction of models are required. The models of one-factor of Schwartz and Heston are defined. The transition density p (x t |x t−1 , θ) becomes intractable, then likelihood-based inference is a difficult problem, a discretized version to the estimation of the parameters and solutions states is written. Geometric Brownian motion models are widely used, due to simplicity and flexibility. Under Geometric Brownian motion (GBM) the price has a lognormal distribution, or equivalently, the log-returns are normally distributed. The GBM for the price at time t of a futures contract with maturity T is denoted F t,T , that is the martingale process, under the risk-neutral probability measure Q where the volatility σ is constant and B t is a Wiener process. Application of Ito's lemma to the no-arbitrage relationship between spot and futures prices, the following representation for the spot price S t is provided where r is the carry cost (including financing, transportation, storage, insurance, etc.) and y is the convenience yield. Also, [11] a relationship between spot and futures prices of a commodity and the convenience yield is established. The authors assume, the convenience yield is a function only of the current spot price, together with the assumption, the risk-free interest rate is a constant ρ. Let S t denote the spot price of a single homogeneous commodity at time t, that follow a stochastic process is where dB t is the increment to a standard Brownian motion, σ t is the instantaneous standard deviation of the spot price, (known) and, µ t is the trend in price, that may be stochastic. The marginal convenience yield C(S t , t) can be defined, that as a function of the current spot price. For tractability, that is assumed C(S t , t) = kS t ; k > 0 (13) simply a proportion of the current spot price, where k is a constant. The performance of the convenience curve can be considered, that is, as an interest rate paid in barrels of oil from borrowing a barrel of oil, then the prices of futures contracts to the delivery of oil can be built. The borrower of a barrel of oil provides storage in the form of inventories to the lender. The lender can be compensated, that for forgoing the benefits associates withholding the barrel of oil. Then, a balance between convenience performance at storage price and shortages of the commodity has existed. Commodity prices by models are determined, convenience returns arise of the interaction between commodity demand, supply, and producer storage decisions. [11] the instantaneous change in the prices of futures contracts is expressed, as Since [25] remark, that in equilibrium only individuals with high convenience yields will hold inventories. The marginal convenience yield, net of storage costs, is inversely proportional to the amount of the commodity held in inventory, in the same way, the marginal convenience yield is directly proportional to the current spot price.
On the other hand, spot price processes through mean-reversion and seasonality can be exhibited, that includes uncertainty about the basis, then [14] the following two-factor process with stochastic convenience yield is introduced where κ is the rate of mean reversion for the convenience yield, α is the mean convenience yield, and λ is the convenience yield risk premium. The fair value relationship between spot and futures under those processes are given, by where Those models for pricing path-dependent options on spot prices are used, that processes display mean-reversion is linked ( see [33], [15]). In particular, two stochastic models to analyze the behavior of commodity prices are studied. In particular, two stochastic models to analyze the behavior of commodity prices are studied, where a novel algorithm to estimate states and parameters is proposed.

The One-Factor Model of Schwartz (1997)
In this section, the one-factor model and derive the corresponding formula for pricing futures contracts are presented, the logarithm of the spot price of a commodity under a mean-reverting process of the Ornstein-Uhlenbeck type is assumed. To develop the one-factor model of [15], the spot price of commodity spot price under a stochastic process is assumed. The one-factor model of [15], the spot price of commodities follow a stochastic process is assumed: To solve the equation analytically (18), introduce the rule of the chain for stochastic differentials is necessitated, the formula of Itô is called, which for x = ln(S), as follow is defined: Applying the formula Itô, the logarithm of the prices by a stochastic process type Ornstein-Uhlenbeck (O-U) can be characterized: where • x t is the net instant convenience performance.
• κ > 0 measures the degree of reversion to the average. • α is the long-term equilibrium level.
• σ characterizes the volatility of the process.
• dB t is a Brownian movement.

ESTIMATION COMMODITY PRICING MODELS
An O-U process under an equivalent measure of martingale is written: • λ is the market risk price.
• dB * t is a Brownian movement under the equivalent measure of martingale. The solution of the equation (23) is given: The conditional distribution (process x in time T ) under the equivalent measure of martingale as a normal is distributed, with mean and variance are specified: • Variance: Note that x T = x(T ) = ln(S(T )) is the spot price of commodities at time T , a random variable log-normal under the equivalent measure of martingale is distributed, with the same parameters. Considering a constant interest rate, the future price of commodities with expiration time T is the expected price of commodities in time T , that an equilibrium measure is maintained, i.e., under an equivalent measure of martingale, consequently: or equivalently: The equation (27) is solution of the partial differential equation: Schwartz's stochastic differential equation (a single factor) in the state-space model is expressed, that completing the partially observed information is allowed; i. e., the unobserved states, and unknown parameters are estimated. The model is expressed: where Kalman filter algorithm for estimating states of the model is usually used, in this article, an alternative algorithm for estimating states and parameters of the model is proposed.

The Heston Model
In this section, the Heston model is introduced. That model by many economic, empirical, and mathematical reasons is chosen. Several important features of low-frequency price data by the model are reproduced: leverage effect, time-varying volatility, fat tails, provides quite reasonable dynamics for the volatility surface, explicit formula for the characteristic function of the asset log-price, consequently, the model calibration procedures are very efficient, and volatility as Brownian diffusion is followed. The Heston model [35] is proposed, the mathematical model to describe the behavior of a bivariate stochastic process is used, a process between the stock price S t and its variance v t is described, that as a generalization of the Black-Scholes option valuation model is initially emerged, assuming the volatility is a stochastic process. The model by the following system of stochastic differential equations is governed: where, the system variables are defined: • S t : is the price of the asset.
• v t : is the volatility of the asset.
• µ : is the expected return on the asset. • θ : is the variation of long-term prices.
• k : is the speed at which volatility tends towards its long-term average.
• ρ : is the correlation of Brownian movements.
• dt = t k − t k−1 : is a small increase in time.
• dB 1,t : is a standard one-dimensional Brownian motion.
• dB 2,t : is a standard one-dimensional Brownian motion.
The discretize the process in the first equation in (31) by Euler discretization is used, the discretization to find the solution of the equation is required. In the following, the process of discretization of s t , ln(s t ) and v t is shown. The stochastic differential equation (SDE) in (31) in integral form can be expressed, as follows Using a discretization of Euler t+dt t µs u du ≈ µs t dt (33) and

ESTIMATION COMMODITY PRICING MODELS
Therefore, the SDE of the stock price is discretized Consider y t = ln (s t ), using Itô's Lemma, Then, Euler's discretization for ln (s t ) is given The discretization of the second equation in (31) v Using a discretization of Euler and where Corr(z u , z s ) = ρ. Finally, the equation of the volatility in the discrete version of Euler is given The Heston's stochastic volatility model as state space models can be written, the discrete version (see [5]): where dt = 1 and y t+1 = lns t+1 are assumed.

Particle Filtering: Algorithm to Estimate States and Parameters
In this section, a sequential Monte Carlo method [34], [2], and recently [23] is proposed, a particle filter, that estimates the states and parameters for long memory stochastic volatility model is based. On smoothing the points of the simulated samples from the parameters and unobserved states of the model is based. The algorithm is consistent and asymptotically normal ( [23]). The procedure in time t is initialized, that generates a combined sample: In this work, a sequential Monte Carlo algorithm [34], [2], and recently [23] is proposed, a particle filter to estimate the states and parameters in long memory stochastic volatility model is based. The particle filter was introduced by [? ] that is estimated the conduct state in nonlinear and non-Gaussian state-space models. The algorithm is defined of the approximating of the filtered density p (x t |y 1:t ) by an empirical distribution that is based on the particles: An equivalent way to define the particle filter is followed: the procedure in time t is initialized, that is generated a combined sample: {x with associated weights: {ω which represents an approximation to the posterior distribution of states and parameters: where δ x0 (x) denotes a mass point of the Dirac delta located at the point x 0 .

Sequential Importance Sampling/Re-sampling Algorithm
Let is considered θ is not fixed; that is, that is varied in time t, for example: To approximate p(θ|y t ), the samples θ (i) t and the weights ω (i) t are generated, the algorithm [29] is used: where m and

ESTIMATION COMMODITY PRICING MODELS
The evolve of states and observations by the equations (30) and (41) is described, in the cases under study in this research. The function of sampling and resampling algorithm of generic importance is summarized: the mass points of the states and parameters are readjusted and states x t , parameters θ t , the associated weights ω t are updated. The joint posterior distribution p(x t , θ t |y 1:t ) is presented. The SISR algorithm has been used to obtain a sample of high-dimensional complex distributions that is occuried in the filtering of non-linear models. The idea of the method is to represent the joint posterior distribution p x t |y 1:t and update θ in conjunction with the x 0:t states.
• Step 1. For t = 1: -For i = 1, . . . , is re-sampled: -The weights are normalized: -Resample: where δ . is the classic Dirac Delta function. • Step 2. For t ≥ 2: and -Thex -For i = 1, . . . , N , is sampled: ω -The weights are normalized: • Step 3. The algorithm outputs: The filtered distribution p(x t |y 1:t ) and the distribution p(θ|y t ), that are approximated: or for and the estimator of θ isθ t :θ the approximation of the combined distribution is given by:

Storvik's SIR Filter
In some situations, the posterior distribution p (θ|x 0:t , y 1:t ) can be depended on a function of the data in terms of the sufficient statistics s t = S (x 0:t , y 1:t ), such that: and the sufficient statistic s t can be recursively updated by:

ESTIMATION COMMODITY PRICING MODELS
The Storvik algorithm is defined by [19] that can be applied to generate the samples of the filtered distribution of the parameters and states. The following decomposition is used, see [20]: p (x 0:t , θ|y 1:t ) ∝ p (x 0:t−1 |y 1:t−1 ) p (θ|s t ) p (x t |x t−1 , θ) p (y t |x t , θ) (71) Storvik's algorithm generates samples from the distribution (71) the following steps are used: with weights w t (i), to obtain a sample from p(x t , s t | y t ).
Steps 2(d)-(e) provide samples from the Bayesian filtering density:

Particle Learning and Smoothing (PLS)
When a particle smoothing algorithm with unknown parameters is used, the states and parameters are estimated that is conditioned to the y 1:t data, that is generated samples x (i) t , θ (i) of the joint posterior distribution p (x t , θ|y 1:t ), where t is the total number of steps. In Carvalho et al. [31] is showed, a filtering algorithm is run, the filtered particles can be resampled to obtain a sample from the smoothed distribution. A decomposition of the smoothed posterior distribution is considered: p (x 0:T , θ|y 1:T ) = p (x T , θ|y 1:T ) The steps of algorithm are: 1. (Forward Filter) Run the particle learning algorithm to generate samples x Step 1 with weights proportional to to generate x

Results
To illustrate the methodology, the states (spot prices) and the parameters of the two stochastic models are estimated, i.e., the one-factor Schwartz model and the Heston model. The predictive capability of the parameter learning algorithm with the models is evaluated. A set of data from the agricultural sector and the mining industry are chosen, the temporary series of contracts with future delivery dates of the primary products soybeans and gold are used. The series consists of 1400 data from August of the year 1990 to September of the year 2017. Two series of data are considered: (1) The data of weekly observations of the prices of contracts with future delivery dates of gold from the New York Mercantile Exchange (NYMEX) is obtained, (2) The data of daily future prices of soybeans from Chicago (CBOT) is obtained, both with expiration dates in the next months, specifically, 1, 6, 12 and 24 months. For both data sets, a total of 4 contracts are considered. The spot price and contracts with a future delivery expiration date are availabled, from (1990). Three sequential Monte Carlo algorithms are Implemented to estimate the states and parameters of the two considered models, the data series of related contracts to the primary products of Soybean and Gold are showed. The first SISR algorithm, that is based on the smoothing of sample points of the latent states and parameters, the procedure can be used to estimate volatility in long-term models dependent on the data. The second algorithm, the states and parameters of the posterior distribution with support in low-dimensional spaces are estimated, that is considered sufficient statistics from the sample. The third algorithm, a learning and smoothing mechanism from the simulated particles with a backward resampling step is used for the choice of weights. The spot price estimate is interest for agents responsible for public policies, central banks, and business decisions in many sectors. The a priori values to the SISR algorithm in the one-factor model of Schwartz are presented: In the Table (1) the estimated posterior means of the parameters by the three algorithms are considered the model of a Schwartz factor is obtained for the adjustment of the Soya series, in addition, a graphic representation in Figure (1) is shown, the converge of the estimated parameters as a function of time is presented. In Figure (2), the unobserved states (spot price) are estimated by the three algorithms, in the green color, the estimate of the posterior mean by the Storvik algorithm are obtained, in the Cyan color, the estimates by the SIRS filter are obtained, in the black color, the estimated value by Carvalho's learning and smoothing particle filter are obtained, and finally the blue color, the real observed price behavior is showed. In Figure (2), the posterior mean values of the estimated 708 ESTIMATION COMMODITY PRICING MODELS latent states by a learning mechanism are obtained that is a correct way because that is the well concentrated values around of the true states.
In Table (2), the estimated parameters of the Heston model for the data series of Soybean are shown. In Figure   Algorithms
(3), a graphic representation of the convergence of the estimated parameters with respect of time is presented. In Figure (4),the unobserved states are estimated by the three proposed algorithms: Storvik (green color), SIRS (Cyan color), Carvalho (black color), and the true states (blue color). In Figure (4) the behavior of the estimated states with respect to the true states are observed, that is the learning structure of the filters with the evolution in time that is followed the same pattern of the true states. The Table (3), the estimates parameters of the one-factor Schwartz model for the Gold series is shown. In Figure   (5), the joint representation of the estimated states by the algorithms are shown: Storvik (green color), SIRS (Cyan color), Carvalho (black color), and the true states (blue color), that can be noted the same similarity with the Soya series, the estimated states are followed the same pattern as the true states.
In Table ( 4), the estimates parameters of the Heston model for the Gold series is reported. In Figure (6 Figure 9. Spot prices and contracts for the one-factor Schwartz model for the Soybean dataset.

Metrics
The performance of the models and the estimation quality of the proposed algorithm is evaluated. The square root of the mean square error (RMSE) as a goodness of fit measure is used.  where • N : denotes the total number of observations. • x i denotes the theoretical spot prices. •x i denotes the estimated spot prices.
In Table (5) The estimated errors by algorithm SISR for one-factor Schwartz model and Heston model are shown, the Soybean and Gold series are considered, the square root of the mean square error (RMSE) measures the quantity of error existing between two data series, in other words, that an estimated value with an observed value is compared. In this article, exist less variability between observed states and estimated states is shown, that validate the results of graphical representations, notably, the algorithm makes a good estimation.

Discussion and Conclusions
In the article, a methodology to estimate the commodity prices from a parameter learning algorithm is proposed. Two stochastic models for estimate spot prices with information of futures contracts are used. The true prices of products in the primary market by the optimal estimation of the models are recovered. The proposed models provide a flexible dynamic, a relationship between observed prices (contracts with a future delivery date) and the unobserved prices (spot price) for data from the agricultural sector and miner is determined. The price volatility of future contracts with maturity dates in the short, medium, and long term by the models in real-time is observed. Seasonal patterns, price peaks, the reversion of the average, price-dependent volatilities, and long-term non-seasonality can be observed; that clear information about the true value of the products is observed. The methodology by real data from the primary sector is validated, the states and parameters simultaneously in realtime are estimated, two commodity price models are presented. A measure the estimation error of the algorithm is estimated, the root of the quadratic mean error is proposed, that result in small estimation errors is observed. Three particle filters to estimate the latent states and parameters are proposed. The proposed methodology is quite general that is applied to a wide class of state space models with non-linear structures and non-Gaussian distributions; in addition, that is allowed the smoothed estimation of states and parameters, the assumptions of linearity and normality of the model are not required, as with the assumptions of the traditional methods (Kalman filter) is used. In this article, the following conditions are proposed, the parameters are not static, in other words, vary over time, and the estimators are obtained with the property of asymptotic normality and are consistent, that is theoretically guaranteed. Future problems with methodology can be solved, more complicated stochastic model structures with non-linear structures and many parameters could be considered. Also, improve the computational efficiency of the algorithm, and consider other variants of the particle filters.