Online Estimation of Time-Varying Volatility Using a Continuous-Discrete LMS Algorithm

The following paper addresses a problem of inference in ﬁnancial engineering, namely, online time-varying volatility estimation. The proposed method is based on an adaptive predictor for the stock price, built from an implicit integration formula. An estimate for the current volatility value which minimizes the mean square prediction error is calculated recursively using an LMS algorithm. The method is then validated on several synthetic examples as well as on real data. Throughout the illustration, the proposed method is compared with both UKF and o ﬄ ine volatility estimation.


INTRODUCTION
In 1973 Black, Scholes and Merton [1,2] reasoned that under certain idealized market assumptions the prices of stocks and the derivatives on these stocks are coupled.One of the crucial assumptions is that the traded asset price S follows dS t = μS t dt + σS t dB t , where B t is a Brownian motion.μ and σ are called, respectively, drift and volatility of the stock; both are deterministic constants.Nevertheless, it turns out that the assumption of constant volatility does not hold in practice.
Traders in the market are supposed to assess returns which have different horizon times in order to predict volatility.Researchers in empirical finance have, therefore, developed an increasing interest in the possibility of uncovering the complex volatility dynamics that exist both within and across different financial markets.Even to the most casual observer of markets, it should be clear that volatility is a random variable.Stochastic volatility models provide a framework for such modeling, especially when dealing with high frequency data.Shephard and Andersen trace the origins of the subject in [3] and attributes it to five sets of people.Back in 1995, the ARCH/GARCH models were a hot topic in econometrics research, and their discoverer, Robert Engle, published a collection of papers on the topic.Now, ten years later, the ARCH/GARCH models are still widely used but their limitations are motivating research into alternative models, specifically, stochastic volatility models (usually abbreviated as SV models).In modern finance, stochastic volatility models represent the latest research which tries to understand financial volatility in continuous time.The resulting process is the nonnegative spot volatility which is assumed to have càdlàg sample paths.The preference given to SV models necessarily follows from the theoretical development of stochastic calculus, which is closely related to continuous time Markov processes.SV models are expected to allow for more comprehensive empirical investigation of the fundamental determinants of certain phenomena: (a) options with different strikes and maturities have different implied volatilities; (b) the empirical distributions of stock returns are leptokurtic.
SV models, consequently, allow for safer measures of risk, for pricing accurately and for hedging options.We refer to Shephard (2005) [4] in order to have a thorough account of the topic of stochastic volatility.All the following studies, for instance, Hull and White (1987) [5], E. M. Stein and J. C. Stein (1991) [6], Heston (1993) [7], Scott (1997), support only offline processing.They aim to calibrate a given model for the volatility dynamics, on the observed sample path of the asset price.The main feature of the method proposed in this paper is an online estimation of volatility: the object to be estimated is one particular trajectory of the volatility process.We use the trajectory of the stock price process, as and when its observation proceeds.Jazwinski in [8] studied the problem of online estimation within continuous time models.In the context of a nonlinear model identification, the use of nonlinear filters such as the unscented Kalman's filter [9,10] is required.
It is proven, however, in [10,11] that traditional UKF is ill-suited for the problem of time-varying volatility estimation.Actually, the UKF never updates prior beliefs, and consequently, it is not able to track volatility fluctuations.We do, however, implement UKF as literature provides no online estimation methods for volatility.Furthermore, we have recourse to an offline estimation method.It is based on an SV model: a continuous time model of volatility dynamics in the form of a stochastic differential equation.Its driving process is Lévy rather than Brownian.The method has been the subject of a recent paper [12].The model frame is built by a "shaping filter" technique [13], using prior information on the covariance function of the squared volatility process.

THE PROPOSED METHOD
To estimate the latent instantaneous volatility σ t of the stock price S t , the stochastic differential equation for the log-price y t = log S t is considered.Applying It ô's formula to (1) yields This SDE may be expressed as The basic idea of the proposed method is to build a predictor from (3) for the observation y t at t = t i+1 .Consequently, (3) is to be discretized at observation instants; this leads us to the question of numerical stability of discretization schemes.It is well known that implicit schemes, such as guarantee numerical stability better [14].Generally, implicit formulae use constant time steps.However, since observations here are made according to arbitrary sampling (i.e., discretization instants are not necessarily equally spaced), only the so-called order-1 and order-2 Adams Moulton formulae are applicable.It is indeed the latter formula (the trapezoidal) that has been chosen: Previously, it has also been used for the identification of a continuous time autoregressive model [15].Equations ( 2)-( 5) lead to The terms holding the Brownian increments ΔB have null expectations.Thus the following predictor y i+1 of the observation y t at t = t i+1 is unbiased: The sense of this choice is that the best model will cause the drift to capture the main course line of the dynamics to the detriment of the diffusion part.Having such a predictor, the estimate of σ i+1 (σ t at t = t i+1 ) that minimizes the mean square prediction error is computed in a recursive way using a stochastic gradient algorithm, the so-called least mean squares algorithm abbreviated to LMS.In this context, the LMS minimizes at each discretization time the following criterion J: using a gradient optimization formula: The resulting formulae are ordered as follows: Initial values y 0 , σ (1)  0 and σ 0 are taken nonstrictly null but arbitrarily small.As usual when using an LMS algorithm, it is the parameter λ that is responsible for the robustness and the right track [16].

ILLUSTRATION
In order to show the performance of the proposed method, different models for the volatility are considered.A constant volatility, for example, is useful in order to evaluate the performance in terms of residual error.A volatility sample path as a step function is interesting in order to evaluate the influence of the initial value on convergence.In addition, it has been widely documented that there is a systematic pattern in average volatility; where this is the case, we will show how estimation of the periodic component of the volatility is feasible.Furthermore, the volatility is modeled as a stochastic process, the solution for an SDE of Vasicek.Finally, we apply our method to real data: the German The proposed method is compared with the UKF for, first, the case of a periodic function of time, and second the case of a "synthetic" stochastic process.UKF is based on a state which has the unobservable volatility process as one of its components.UKF equations of the time and measurement updates for the first moment μ of the conditional density are, respectively, Estanding for the mathematical expectation.UKF, thus, does not update prior estimates μ(t i+1 |t i ), and consequently it is not able to track time-varying volatility.Similar behavior is exhibited in [10,11].
Next, a comparison is made between the above method and an offline estimation of the volatility.The latter was proposed in [12] which deals with the construction of a black-box continuous time model for the squared volatility process in the form of a stochastic differential equation.The starting point in this construction is a parametric form for the covariance function of this process.The model frame derives from a control theory technique known as the shaping filter.We give a brief account of the work presented in [12] and show that our present study outperforms it.
As regards observations, they are made according to both periodic and nonperiodic sampling schemes.For instance, the case of jitter sampling, as in [15], is considered in Section 3.2.The obtained performance is as good as that of a periodic sampling scheme.

Constant volatility
The observations are simulated with a volatility of 0.15.
The initial value of the volatility, in the proposed method, is deliberately taken equal to the true value (= 0.15) so that we evaluate the residual estimation error.A periodic sampling scheme has been used.The result is reported in Figure 1.Both the mean value and the standard deviation of the relative error of estimation are about 1% and 6%, respectively.They are calculated by time averaging since the volatility value is constant along its trajectory.

Volatility as a step function
In order to illustrate the convergence behavior of the proposed method, a step function with the initial value of 0.1 and the final value of 0.2 is taken as the volatility sample path.The proposed method is implemented with an initial value of 0.1 for the volatility.A jitter sampling scheme has been used with maximum value of half the sampling period.Many simulations have been carried out with different values of λ; the value 0.04 for λ makes a good tradeoff between robustness and right track.The result is reported in Figure 2; it shows the capability of the algorithm to follow rapid variations even for nonuniformly sampled data.Both the mean value and the standard deviation of the relative error of estimation are about 1% and 10%, respectively.Here again they are calculated by time averaging; this is legitimate since there is piecewise repetition of the volatility value along its trajectory.To explore further the performance evaluation of this result, we have computed the Theil index.It is approximately 3  Here N is the number of samples in the reference trajectory to be estimated.σ est is the estimate of the volatility σ t at t = t i+1 , denoted by σ i+1 in Section 2, and σ ref is the reference: the (true) volatility σ t at t = t i+1 , denoted σ i+1 (i = 0, . . ., N − 1).In addition, Monte Carlo simulations have been carried out: the mean sample path for 100 estimated trajectories of the volatility is reported in Figure 3.The mean value and the standard deviation of its relative error of estimation are about 1% and 5%, respectively.This shows that the standard deviation of the estimation error drops significantly as the simulation number increases.That is, as expected, the empirical mean sample paths are to converge to the true mean.
As has been said in the introduction to Section 3, the parameters can be chosen arbitrarily within all synthetic examples.The only essential thing to take into consideration is the realistic values of the volatility.The general validity of our method should thus be studied by varying these parameters.They are the initial and the final values of the step function in the context of this subsection.Column 1 in Table 1 shows initial values of three different step functions; column 2 shows their corresponding final values.Columns 3 and 4 show the mean value and the standard deviation of the relative error of estimation obtained by Monte Carlo simulations (25 estimated trajectories of the volatility for each couple of parameters).The last two columns show the mean Theil index of the 25 estimated trajectories of the volatility using our method versus the Theil index of UKF for each couple of parameters.

Volatility as a deterministic periodic function of time
Whenever the volatility is subject to seasonality, we wish to recover the season(s) using our method.We consider the following deterministic function of time for the volatility  trajectory: The pulsations ω 1 and ω 2 correspond to a one-week and a one-day seasonality; this is, for instance, the case of German price treated in 3.6.a 0 , a 1 and a 2 are chosen so as to have realistic values of the volatility.In the simulation of Figure 4, they are 0.15, 0.05, and 0.01, respectively.Both the true volatility and its estimate for a periodic sampling scheme, and for λ of 0.07, are plotted in Figure 4.The estimated volatility using UKF is constant, yet the proposed method is able to track the volatility oscillations.The Theil index is about 10 −3 ; UKF yields a Theil index of 10 −2 .The mean value and the standard deviation of the relative error are about 1% and 16%, respectively; they are calculated by time averaging.To justify this, we do check error ergodicity.This is done by fixing an instant, repeating again the simulation several times with respect to the same volatility trajectory till this instant.The mean value and the standard deviation of the relative error for this instant are obtained by averaging on simulations.Their values are in the order of what is given above.Besides, we proceed likewise in the following (the mean value and the standard deviation of the relative error are to be calculated by time averaging).
The mean trajectory of 100 estimated trajectories of the volatility is reported in Figure 5.The mean value and the standard deviation of its relative error of estimation are about 1% and 8%, respectively.In addition, the power spectral densities (PSD) for the true volatility sample path and the mean of its estimates are confronted in Figure 6; the two PSDs therein are clearly close to each other.
We furthermore vary the parameters a 1 and a 2 and perform Monte Carlo simulations (100 estimated trajectories of the volatility for each couple of parameters) so that we obtain the results in Table 2.

Volatility as a stochastic process
To synthesize sample paths of the volatility process as well as the stock price, the following SDE of Vasicek is considered: where α = 0.0001, θ = 0.15, and ξ = 0.0007.We assume the drift μ is known (μ = 0.015).The true volatility sample path and the estimated one, using both the proposed method and UKF, are reported in Figure 7.The volatility is estimated at every half hour for 416 days.For this simulation, we choose the initial value of the volatility equal to θ(= 0.15).As above, the estimated volatility using UKF is constant.The proposed method, however, is able to track the volatility fluctuations.The empirical distribution of the estimation error for the sample path in Figure 7 is reported in Figure 8. Like UKF, the proposed method is subject to bias, but the bias is clearly smaller.The standard deviation obtained with UKF is 0.033, whereas within the proposed method, it is 0.015.path exhibits a volatility clustering phenomenon: periods of high-price fluctuations are followed by periods of high fluctuations, and the same can be said about periods of lowprice fluctuations.The implementation result on this sample path is shown in Figure 10.Notice the beginning of a period of high volatility around the 700th day; this corresponds to the Asian financial crisis of October 1997.

Comparison with offline estimation of the volatility
We assume prior information about the unknown process (σ t ) 2 : its stationarity in the large sense and a parametric model for its covariance function.Let the covariance function of the process (σ t ) 2 be given by the following formula: where D is the process variance.This type of covariance function allows one to fit the observed time dependence in the returns.Such a covariance function includes memory in the correlation pattern of the volatility.The spectral density of (σ t ) 2 is then given by the formula   The spectral density s(ω) is rewritten as where represents the transfer function of a stationary linear system; the system is, furthermore, stable as the root of F(s) is in the left half-plane of the complex variable s.Recalling that 1/2π is the spectral density of a white noise of intensity 1, we  come to the following conclusion.(σ t ) 2 may be considered as the response of the filter whose transfer function is Φ(s) to a white noise with unit intensity.From the ordinary differential equation describing such a filter, we obtain the following stochastic differential equation as a model for the squared volatility process (σ t ) 2 .This is the first state component denoted by X 1 t : (20) Here, W is a stochastic process with independent and stationary increments of intensity 1.If we suppose that W starts at 0 and that its trajectories are continuous in probability, then we can give it the name Lévy process.We suppose further the existence of stationary solutions to the SDE when W has positive increments so as to assure the

Figure 3 :
Figure 3: True volatility (dashed) versus the mean for 100 of its estimates (continuous).

Figure 5 :
Figure 5: True volatility (dashed) versus the mean for 100 of its estimates (continuous).

Figure 9 Figure 6 :
Figure 9 shows the daily price of the Hang Seng index of the Hong Kong market from 1995 to 2007.This sample

Figure 8 :
Figure 8: Empirical distribution for the estimation error.(a) The proposed method, (b) UKF.