Nonlinear Time Series Modeling:
A Unified Perspective, Algorithm, and Application
Subhadeep Mukhopadhyay
Department of Statistical Science, Temple University
Emanuel Parzen
Department of Statistics, Texas A&M University
Abstract

arXiv:1308.0642v4 [] 24 Dec 2017

A new comprehensive approach to nonlinear time series analysis and modeling is
developed in the present paper. We introduce novel data-specific mid-distribution
based Legendre Polynomial (LP) like nonlinear transformations of the original time
series {Y (t)} that enables us to adapt all the existing stationary linear Gaussian time
series modeling strategy and made it applicable for non-Gaussian and nonlinear processes in a robust fashion. The emphasis of the present paper is on empirical time series
modeling via the algorithm LPTime. We demonstrate the effectiveness of our theoretical framework using daily S&P 500 return data between Jan/2/1963 - Dec/31/2009.
Our proposed LPTime algorithm systematically discovers all the ‘stylized facts’ of the
financial time series automatically all at once, which were previously noted by many
researchers one at a time.

Keywords and phrases: Nonparametric time series modeling, Nonlinearity, Unified time
series algorithm, Exploratory diagnostics.

Contents
1 Introduction

2

2 From Linear to Nonlinear Modeling

3

3 Nonparametric LPTime Analysis
3.1 The Data and LP-Transformation . . . . . . . . . . . . . . . .
3.2 Marginal Modeling . . . . . . . . . . . . . . . . . . . . . . . .
3.2.1 Non-normality Diagnosis . . . . . . . . . . . . . . . . .
3.3 Copula Dependence Modeling . . . . . . . . . . . . . . . . . .
3.3.1 Nonparametric Serial Copula . . . . . . . . . . . . . .
3.3.2 LP-Comomemnt of lag h . . . . . . . . . . . . . . . . .
3.3.3 LP-Correlogram, Evidence and Source of Nonlinearity
3.3.4 AutoLPinfor: Nonlinear Correlation Measure . . . . .
3.3.5 Nonparametric Estimation of Blomqvist’s Beta . . . .
3.3.6 Nonstationarity Diagnosis, LP-Comoment Approach .
3.4 Local Dependence Modeling . . . . . . . . . . . . . . . . . . .
3.4.1 Quantile Correlation Plot and Test for Asymmetry . .
3.4.2 Conditional LPinfor Dependence Measure . . . . . . .

1

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.

4
5
5
5
6
6
8
8
9
10
11
12
12
12

3.5
3.6
3.7

Non-Crossing Conditional Quantile Modeling . . . . . . . . . . . . . . . . . . . . .
Nonlinear Spectrum Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Nonparametric Model Specification . . . . . . . . . . . . . . . . . . . . . . . . . . .

14
16
18

4 Conclusion

18

5 References

20

1

Introduction

When one observes a sample Y (t), t = 1, . . . , T , of a (discrete parameter) time series Y (t),
one seeks to nonparametrically learn from the data a stochastic model with two purposes:
(a1) scientific understanding; (a2) forecasting (predict future values of the time series under
the assumption that the future obeys the same laws as the past). Our prime focus in this
paper is on developing a nonparametric empirical modeling technique for nonlinear (stationary) time series that can be used by data scientists as a practical tool for obtaining insights
into (i) the temporal dynamic patterns and (ii) the internal data generating mechanism–a
crucial step for achieving (a1) and (a2).
Under the assumption that the time series is stationary (which can be extended to asymptotically stationary) the distribution of Y (t) is identical for all t, the joint distribution of
Y (t) and Y (t + h) depends only on lag h. Typical estimation goals are as follows:
(1) Marginal modeling. The identification of marginal probability law (in particular, the
heavy tailed marginal densities) of a time series plays a vital role in financial econometrics.
Notations: Common quantile Q, inverse of distribution function F , respectively denoted
Q(u; Y ), 0 < u < 1 and F (y; Y ). Mid-distribution is defined as F mid (y; Y ) = F (y; Y ) −
.5 Pr(Y (t) = y).
(2) Correlation modeling. Covariance function (defined for positive and negative lag h)
R(h; Y ) = Cov[Y (t), Y (t+h)]. R(0; Y ) = Var[Y (t)], µ = E[Y (t)] assumed 0 in our prediction
theory. Correlation function ρ(h) = Cor[Y (t), Y (t + h)] = R(h; Y )/R(0; Y ).
(3) Frequency-domain modeling. When covariance is absolutely summable, define spectral
P
density function f (ω; Y ) = R(h; Y ) e−2πiωh , −1/2 < ω < 1/2.
(4) Time-domain modeling. Time domain model is linear filter relating Y (t) to white noise
(t), N (0, 1) independent random variables. Autoregressive scheme of order m, a predominant linear time series technique for modeling conditional mean is defined as (assuming
E[Y (t)] = 0 )
Y (t) − a(1; m)Y (t − 1) − . . . − a(m; m)Y (t − m) = σm (t),

2

(1.1)

with the spectral density function given by
f (ω; Y ) =

2
σm

1−

Pm

k=1

a(k; m)e2πiωk

2.

(1.2)

To fit an AR model, compute linear predictor of Y (t) given Y (t − j), j = 1, . . . , m by


Y µ,m [t] = E Y (t) | Y (t−1), . . . , Y (t−m) = a(1; m)Y (t−1)+· · ·+a(m; m)Y (t−m). (1.3)
Verify that the prediction error Y [t] − Y µ,m [t] are white noise. Best fitting AR order is
identified by Akaike criterion AIC (or Schwarz’s criterion BIC) as value of m minimizes
AIC(m) = 2 log σm + 2m/n.
In what follows we aim to develop a parallel modeling framework for non-linear time series.

2

From Linear to Nonlinear Modeling

Our approach to nonlinear modeling, called LPTime, is via approximate calculation of conditional expectation E[Y (t)|Y (t − 1), . . . , Y (t − m)]. Because with probability one Q(F (Y )) =
Y , one can prove that the conditional expectation of Y (t) given past values Y (t − j) is equal
to (with probability one) conditional expectation of Y (t) given past values F mid (Y (t −
j)), which can be approximated by linear orthogonal series expansion in score functions
Tk [F mid (Y (t − j))] constructed by Gram Schmidt orthonormalization of powers of
T1 =

F mid (Y (t); Y ) − .5
,
σ[F mid (Y (t); Y )]

(2.1)

where σ[F mid (Y (t);q
Y )] si the standard deviation of the mid-distribution transform random
P
variable given by (1 − y p3 (y))/12, p(y) denotes the probability mass function of Y .
This score polynomials allows us to simultaneously tackle the discrete (say count-valued)
and continuous time series. Note that for Y is continuous, The T1 reduces to:
T1 =

√

12 F (Y (t)) − .5 .

(2.2)

and all the higher-order polynomials Tj can be compactly expressed as Legj [F (Y )], where
Legj (u), 0 < u < 1 denotes orthonormal Legendre polynomials. It is worthwhile to note
that Tj are orthonormal polynomials of mid-ranks (instead of polynomials of the original
y’s), which injects robustness into our analysis while allowing to capture nonlinear patterns.
Having constructed score functions of y denoted by Tj , we transform into unit interval by
letting y = Q(u; Y ) and defining




Sj (u; Y ) = Tj F mid (Q(u; Y )) , Tj (y; Y ) = Sj F mid (Y (t)) .

(2.3)

In general our score functions are custom constructed for each distribution function F which
can be discrete or continuous.
3

−10
0.0 1.0 −20
0.5 1.5
0 1 2 −1.0
1
−1

0

S4

2

3 −2

S3

S2

−1.5

S1

Z

0

10

LP−Transformed Time Series

Time

Figure 1: LP-transformed S&P 500 daily stock returns between Oct 1986 and Oct 1988.
This is just a small part of the full time series from Jan/2/1963 - Dec/31/2009 (cf. Sec 3.1).

3

Nonparametric LPTime Analysis

Our LPTime empirical time series modeling strategy to nonlinear modeling of a univariate
time series Y (t) is based on linear modeling of the multivariate time series:
Vec(YS)(t) =



T
YS1 (t), . . . , YSk (t) ,

(3.1)

where YSk (t) = Tk [F mid (Y (t))], our tailor-made orthonormal mid-rank based nonlinear
transformed series. We summarize below the main steps of the algorithm LPTime. To
better understand the functionality and applicability of LPTime we break it into several
inter-connected steps each of which highlights:
(a) Algorithmic modeling aspect [how it works]
(b) Required theoretical ideas and notions [why it works]
(c) Application to daily S&P 500 return data between Jan/2/1963 - Dec/31/2009 [empirical
proof-of-work].
4

1.5

1.5

0.5

1.0

0.4

1.0
0.5

0.5
−15

−10

−5

0

5

10

0.0

0.0

Density

0.3
0.2
0.1
0.0

−20

0.0

(a)

0.2

0.4

0.6

0.8

1.0

0.0

(b)

0.2

0.4

0.6

0.8

1.0

(c)

Figure 2: (a) The marginal distribution of daily returns; (b) plots the histogram of Φ(yi ) and
display the LP-estimated comparison density curve. and (c) shows the associated comparison
density estimate with G as t-distribution with 2 degrees of freedom.

3.1

The Data and LP-Transformation

The data used in this paper is daily S&P 500 return data between Jan/2/1963 - Dec/31/2009
(defined as log(Pt /Pt−1 ) where Pt is the closing price on trading day t). We begin our
modeling process by transforming the given univariate time series {Y (t)} into multiple
(robust) time series by means of a special data-analytic construction rules described in
(2.1)-(2.3) and (3.1). We display the original “normalized” time series Z(Y (t)) = (Y (t) −
E[Y (t)])/σ[Y (t)] and the transformed time series YS1 (t), . . . , YSk (t) on a single plot.
Fig 1 shows the first look at the transformed S&P 500 return data between Oct 1986 and Oct
1988. These newly constructed time series works as a universal preprocessor for any time
series modeling in contrast with other adhoc power transformations. In the next sections
we will describe how the temporal patterns of these multivariate LP-transformed series
Vec(YS)(t) = {YS1 (t), . . . , YSk (t)} generate various insights for the time series {Y (t)} in a
organized fashion.

3.2

Marginal Modeling

Our time series modeling starts with the identification of probability distributions.
3.2.1

Non-normality Diagnosis

Does the Normal probability distribution provide a good fit to the S&P 500 return data?
Fig 2(a) clearly indicates distribution of daily return is certainly non-normal. At this point
the natural question is how the distribution is different from the assumed normal? A quick
insight into this question can be gained by looking at the distribution of the random variable
U = G(Y ), called comparison density (Parzen, 1997, Mukhopadhyay, 2017) given by:
d(u; G, F ) =

f (Q(u; G))
,
g(Q(u; G))
5

0 ≤ u ≤ 1,

(3.2)

where Q(u; G) = inf{x : G(x) ≥ u} is the quantile function. The flat uniform shape of
the estimated comparison density provides a quick graphical diagnostic to test the fit of
the parametric G to the true unknown distribution F . The Legendre polynomial based
orthogonal series comparison density estimator is given by
X
d(u; G, F ) = 1 +
LP[j; G, F ] Legj (u), 0 < u < 1

(3.3)

j

where the Fourier coefficients LP[j; G, F ] = E[Legj ◦ G(Y )].
For G = Φ, Fig 2(b) displays the histogram of Ui = Φ(Yi ) for i = 1, . . . , n. The corresponding
b G, F ) = 1 − 0.271 Leg (u) − 0.021 Leg (u) + 0.193 Leg (u)
comparison density estimate d(u;
2

3

4

is shown in blue curve, which reflects the fact that distribution of daily return has (i)
sharp peaked (inverted “U” shape), (ii) negatively skewed with (iii) fatter tails than the
Gaussian distribution. We can carry our similar analysis by asking whether t-distribution
with degrees of freedom 2 provides a better fit. Fig 2(c) demonstrates the full analysis,
b G, F ) = 1 − 0.492 Leg2 (u) − 0.015 Leg3 (u) +
where the estimated comparison density d(u;
0.084 Leg4 (u) indicates (iv) t-distribution fits the data better than normal, specially in the
tails, although not a fully adequate model.
The shape of comparison density (along with the histogram of Ui = G(Yi ), i = 1, . . . , n)
captures and exposes the adequacy of the assumed model G for the true unknown F –thus
act as an exploratory as well as confirmatory tool.

3.3

Copula Dependence Modeling

Distinguishing uncorrelatedness and independence by properly quantifying association is an
essential task in empirical nonlinear time series modeling.
3.3.1

Nonparametric Serial Copula

We display the nonparametrically estimated smooth serial copula density cop(u, v; Y (t), Y (t+
h)) to get much finer understanding of the lagged interdependence structure of a stationary
time series. For continuous distribution define the copula density for the pair (Y (t), Y (t+h))
as the joint density of U = F (Y (t)) and V = F (Y (t + h)), which is estimated by sample
mid-distribution transform Ũ = Femid (Y (t)), Ṽ = Femid (Y (t + h)). Following Mukhopadhyay and Parzen (2014) and Parzen and Mukhopadhyay (2012), we expand copula density
(square integrable) in a orthogonal series of product LP-basis functions as
X
cop(u, v; Y (t), Y (t+h)) − 1 =
LP[j, k; Y (t), Y (t+h)] Sj (u; Y (t)) Sk (v; Y (t+h)), (3.4)
j,k

where Sj (u; Y (t)) = YSj (Q(u; Y (t)); Y (t)). The Eq 3.4 allows us to pictorially represent the
information present in the LP-comoment matrix via copula density. The various “shapes”
of copula density gives insight into the structure and dynamics of the time series.
6

Now we apply this nonparametric copula estimation theory to model the temporal dependence structure of S&P return data. The copula density estimate cd
op(u, v; Y (t), Y (t + 1))
based on the smooth LP comoments is displayed in Fig 3. The shape of the copula density shows strong evidence of asymmetric tail-dependence. Note that the dependence is only
present in the extreme quantiles - another well-known stylized fact of economic and financial
time series.

3
2
1

0.034

1.0
0.8
0.6
0.4
0.2
0.00.0

●

0.4

0.2

●

0.6

1.0

0.8

●

0.032

●
●

0.030

●

●

●

0.028

●

0.026

●

●

0.024

●
●

0.022

●
●
●

1

3

6

9

12

15

Figure 3: Top: Nonparametric smooth serial copula density (lag 1) estimate of S&P return
data. Bottom: BIC plot to select the significant LP-comoments computed in Eq. 3.5.

7

3.3.2

LP-Comomemnt of lag h

Here we will introduce the concept of LP-comoment to get a complete understanding of the
nature of serial dependence present in the data. LP-comoment of lag h is defined as the
joint covariance of Vec(YS)(t) and Vec(YS)(t + h).
The lag 1 LP-comoment matrix for S&P 500 return data is

0.0705∗ −0.0617∗

 0.0074


0.1542∗
LP Y (t), Y (t + 1) = 
−0.0104 −0.0071

0.0166
0.0438∗

displayed below,
0.0199

0.0113




0.0652∗ 

0.0262 −0.0355

∗
0.0113 0.0698
0.0077

(3.5)

To identify the significant elements we first rank order the squared LP comoments. Then we
take the penalized cumulative sum of m comoments using BIC criterion 2m log(n)/n, n is
sample size and choose the m for which BIC is maximum. The complete BIC path for S&P
500 data is shown in figure 3, which selects top six comments also denoted by ∗ in the LPcomoment matrix display (Eq. 3.5). By making all those uninteresting “small” comoments
c . The linear
equal to zero we get the “smooth” LP Comoment matrix denoted by LP
auto-correlation is capture by the LP[1, 1; Y (t), Y (t + 1)] = E[YS1 (t) YS1 (t + 1)] term. The
presence of higher order significant terms in the LP comoment matrix indicate the possible
nonlinearity. Another interesting point to note that CORR[Y (t), Y (t + 1)] = .027, where as
the auto-correlation between the mid-rank transformed data CORR[F mid (Y (t)), F mid (Y (t +
1))] = .071, considerably larger and picked by BIC criterion. This is an interesting fact as
it indicates rank-transform time series (YS1 (t)) is much more predictable than the original
raw time series Y (t).
3.3.3

LP-Correlogram, Evidence and Source of Nonlinearity

We provide nonparametric exploratory test for (non)linearity (spectral domain test is given
in Sec 3.6). Plot the correlogram of YS1 (t), . . . , YS4 (t): (a) diagnose possible nonlinearity,
and (b) identify possible sources. This constitutes the important building block for methods
of model identification. LP-correlogram generalizes the classical sample autocorrelation
function (ACF). Applying the acf() R function on Vec(YS)(t) generates the graphical
display of our proposed LP-correlogram plot.
Fig 4 shows the LP-correlogram of S&P stock return data. Panel A shows the absence of
linear autocorrelation which is known as a efficient market hypothesis in finance literature.
A prominent auto-correlation pattern for the series YS2 (t) (panel B Fig 4) is the source
of nonlinearity. This fact is known as “volatility clustering”, which says that large price
fluctuation is more likely to be followed by large price fluctuations. Also the slow decay of
the autocorrelation of the series YS2 (t) can be interpreted as an indication of long-memory
volatility structure.
8

0.0

0.0

0.2

0.2

0.4

0.4

0.6

0.6

0.8

0.8

1.0

S2

1.0

S1

50

100

150

0

50
Lag

S3

S4

100

150

100

150

0.8
0.6
0.4
0.2
0.0

0.0

0.2

0.4

0.6

0.8

1.0

Lag

1.0

0

0

50

100

150

0

50

Lag

Lag

Figure 4: LP-correlogram: sample autocorrelations of LP-Transformed time series. The
decay rate of the sample autocorrelations of YS2 (t) appears to be much slower than the
exponential decay of ARMA process, implying possible long memory behaviour.
3.3.4

AutoLPinfor: Nonlinear Correlation Measure

We display the sample AutoLPinfor plot - diagnostic tool for nonlinear autocorrelation. We
define the lag h AutoLPinfor as the squared Frobenius norm of the smooth-LP comoment
matrix of lag h,
AutoLPinfor(h) =

2

X

LP[j, k; Y(t), Y(t + h)] ,

(3.6)

j,k

where sum is over BIC selected j, k for which LP comoments are significantly non-zero.
Our robust nonparametric measure can be viewed as capturing the deviation of copula
density from uniformity:
ZZ
AutoLPinfor(h) =

cop2 [u, v; Y(t), Y(t + h)] du dv − 1,

(3.7)

which is closely related to the entropy measure of association proposed in Granger and Lin
(1994)
ZZ
Granger–Lin(h) =

cop[u, v; Y(t), Y(t + h)] log cop[u, v; Y(t), Y(t + h)] du dv.
9

(3.8)

ACF Plot

0.0

0.02

0.2

0.04

0.4

0.06

0.6

0.8

0.08

1.0

AutoLPinfor: Nonlinear Correlation Detector

0

50

100

150

0

Lag in days

50

100

150

Lag in days

Figure 5: Left: ACF plot of S&P 500 data; Right: AutoLPinfor Plot upto lag 150.
It can be showed using Taylor series expansion that asymptotically
AutoLPinfor(h) ≈ 2 × Granger–Lin(h).

(3.9)

An excellent discussion on the role of information theory methods for unified time series
analysis is given in Parzen (1992) and Brillinger (2004). For an extensive survey of tests of
independence for nonlinear processes see Chapter 7.7 of Terasvirta et al. (2010). AutoLPinfor is a new information theoretic nonlinear autocorrelation measure which detects generic
association and serial dependence present in a time-series. Contrast the AutoLPinfor plot
for S&P 500 return data shown in Fig 5 with the acf plot (left panel). This underlies the
need for building a nonlinear time series model, which we will be discussing next.
3.3.5

Nonparametric Estimation of Blomqvist’s Beta

Estimate the Blomqvist’s β (also known as medial correlation coefficient) of lag h by using
the LP-copula estimate in the following equation,
Z1/2Z1/2
βbLP (h; Y (t)) := − 1 + 4
0



cd
op u, v; Y (t), Y (t + h) du dv

(3.10)

0

The β values −1, 0 and 1 interpreted as reverse correlation, independence and perfect correlation, respectively. Note that,
Blomqvist’s β
AutoLPinfor

: Normalized distance of copula distribution Cop(u, v) from independence copula uv
: Distance of copula density cop(u, v) from uniformity 1.

For S&P 500 return data we compute the following dependence numbers,
10

3
2
1

1.0
0.8
0.6
Y[
t] 0.4
0.2
0.00.0

0.4

0.2

0.6

1.0

0.8

Time

Figure 6: LP copula diagnostic for detecting non-stationarity in S&P 500 return data.
βbLP (1; Y (t)) = 0.0528
βbLP (1; YS1 (t)) = 0.0528
βbLP (1; YS2 (t)) = 0.0729
βbLP (1; YS3 (t)) = 0.0
βbLP (1; YS4 (t)) = 0.003.
3.3.6

Nonstationarity Diagnosis, LP-Comoment Approach

Viewing the time index T = 1, . . . , n as covariate we propose a nonstationarity diagnosis
based on LP-comoments of Y (t) and the time index variable T . Our treatment has the
ability to detect time varying nature of mean, variance, skewness and so on represented by
various custom-made LP-transformed time series.
For S&P data we computed the following LP comoment matrix to investigate the nonstationarity:


0.180∗

0.012

−0.010

0.058∗





−0.005 −0.034 −0.036 0.080∗ 




LP T, Y (t) = 

∗
−0.016
0.115
0.001
−0.001


∗
0.024 −0.040 −0.010 0.049

(3.11)

This indicates presence of slight non-stationarity behavior of variance or volatility (Y S2 (t))
and kurtosis of tail-thickness(Y S4 (t)). Similar to AutoLPinfor we propose the following
statistic for detecting nonstationarity
LPinfor[Y(t), T] =

X
j,k

11

2

LP[j, k; T, Y(t)] .

(3.12)

We can also generate the corresponding smooth copula density of (T, Y (t)) based on smooth


LP T, Y (t) matrix to visualize the time varying information as in Fig 6.

3.4
3.4.1

Local Dependence Modeling
Quantile Correlation Plot and Test for Asymmetry

Display the quantile correlation plot, a copula distribution based graphical diagnostic to
visually examine the asymmetry of dependence. The goal is to get more insight into the
nature of tail-correlation.
Motivated by the concept of lower and upper tail-dependence coefficient we define the Quantile Correlation Function (QCF) as the following in terms of copula distribution function of
(Y (t), Y (t + h)) denoted by Cop(u, v; Y (t), Y (t + h)) := Cop(u, v; h),


Cop(u, u; h) 
I
λ u; Y (t), Y (t + h) :=
u≤.5
u

+

1 − 2u + Cop(u, u; h) 
I
,
u>.5
1−u

(3.13)

Our nonparametric estimate of quantile correlation function is based on the LP-copula den

bLP u; Y (t), Y (t + h) . Fig 7 shows the corresponding quantile
sity which we denote as λ
correlation plot for S&P 500 data. The dotted-line represent QCF under Independence assumption. Deviation from this line helps us to better understand the nature of asymmetry.
bG [u; Y (t), Y (t + h)] using the fitted Gaussian copula
We compute λ
d G (u, v; Y (t), Y (t + h)) = Φ(Φ−1 (u), Φ−1 (v); Σ̂ = S)
Cop

(3.14)

where S is sample covariance matrix. The dark-green line in Fig 7 shows the corresponding
curve which is almost identical with the “no dependence” curve, albeit misleading. The
reason is Gaussian copula is characterized by linear correlation, while S&P data is highly
nonlinear in nature. As the linear auto-correlation of stock return is almost zero, we have
approximately Φ(Φ−1 (u), Φ−1 (u); Σ̂ = S) ≈ Φ(Φ−1 (u))Φ(Φ−1 (u)) = u2 . Similar to Gaussian
copula there are several others parametric copula families, which can give similar misleading
conclusions. This simple illustration reminds us the pernicious effect of not “looking into
the data”.
3.4.2

Conditional LPinfor Dependence Measure

For more transparent and clear insights into the asymmetric nature of the tail dependence
we need to introduce the concept of conditional dependence. In what follows we propose a
conditional LPinfor function LPinfor(Y (t+h)|Y (t) = Q(u; Y (t))) - quantile based diagnostic
for tracking how the dependence of Y (t + h) on Y (t) changing at various quantiles.
To quantify the conditional dependence, we seek to estimate f (y; Y (t + h)|Y (t))/f (y; Y (t +
h)). A brute force approach estimates separately the conditional distribution and the unconditional distribution and take the ratio to estimate this arbitrary function. An alternative
12

0.0

0.1

0.2

0.3

0.4

0.5

LP−estimate
Gaussian copula
Independence

0.0

0.2

0.4

0.6

0.8

1.0

Quantiles

bLP [u; Y (t), Y (t + 1)]. Detects
Figure 7: Estimated Quantile Correlation Function (QCF) λ
asymmetry in the tail-dependence between lower-left quadrant and upper-right quadrant for
S&P 500 return data. The red dotted line denotes the quantile correlation function under
dependence. The dark-green line shows the quantile correlation curve for fitted Gaussian
copula.
elegant way is to recognize that by “going to quantile domain” (i.e., Y (t+h) = Q(v; Y (t+v))
and Y (t) = Q(u; Y (t))) we can interpret the ratio as “slices” of copula density, which we
call conditional comparison density:
X


d v; Y (t + h), Y (t + h)|Y (t) = Q(u; Y (t)) = 1 +
LP[j; h, u] Sj (v; Y (t + h)),

(3.15)

j

where the LP-Fourier orthogonal coefficients LP[j; h] are given by
X
LP[j; h, u] =
LP[j, k; Y (t), Y (t + h)] Sk (u; Y (t)).
k

Define conditional LPinfor as
X


2
LPinfor Y (t + h)|Y (t) = Q(u; Y (t)) =
LP[j; h, u] .

(3.16)

j

We use this theory to investigate the conditional dependency structure of S&P 500 return
data. Fig 8(a) traces out the complete path of the estimated LPinfor[Y (t + h) | Y (t) =
Q(u; Y (t)] function, which indicates the high asymmetric tail-correlation. This conditional
correlation curves can be viewed as “local” dependence measure. An excellent discussion on
this topic is given in Section 3.3.8 of Terasvirta, Tjostheim, and Granger (2010).
13

●

0.0

●

0.0

●●●

●●●●●

0.2

●●●●●●●

●●

●●●

●●●●●

0.4

●●●

0.6

Quantiles of Y(t)

●●●●●●●●●●●●●●

0.8

●

1.0

●

0.0

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

0.6

●

●

●

0.4

●

●

●

●
●
●

0.2

●

●

●

●

●

Volatility level

●

●

●

●

●

●

●
●

●
●

●

●●

●●

●

●
●

●

●

●
●

●

●
●
●
●
●

−0.2

●

−0.10

0.1

●

●

●

●

●

●

●

0.0

Mean level
●

−0.05

0.3

●

0.00

0.05

0.5
0.4

●

0.2

Conditional LPinfor

0.6

0.10

●

0.2

0.4

0.6

0.8

1.0

0.0

0.2

Quantiles of Y(t)

0.4

●

●●

●●●

●

●

●

●

●

0.6

●

●

●●

●●●●●●●

0.8

1.0

Quantiles of Y(t)

Figure 8: (a) The conditional LPinfor curve is shown for the pair [Y (t), Y (t + 1)]. Clearly
shows the asymmetric dependence in the tails and almost nothing is going on in between. (b) and (c) display how the mean and volatility levels of conditional distribution
f [y; Y (t + 1)|Y (t) = Q(u; Y (t))] changes with respect to the unconditional marginal distribution f (y; Y (t)) at different quantiles.
At this point we can legitimately ask what aspects of the conditional distributions are changing most? Fig 8 (b,c) display only the two coefficients LP[1; h, u] and LP[2; h, u] for the S&P
500 return data for the pairs (Y (t), Y (t + 1)). These two coefficients represent how the mean
and the volatility levels of the conditional density changing with respect to the unconditional
reference distribution. The typical asymmetric shape of conditional volatility shown in right
panel of Fig 8 (b,c) indicates what is known as the “leverage effect” - future stock volatility
negatively correlated with past stock return, i.e. stock volatility tends to increase when
stock prices drop.

3.5

Non-Crossing Conditional Quantile Modeling

Display the nonparametrically estimated conditional quantile curves of Y (t + h) given Y (t).
b h, u) to
Our new modeling approach uses the estimated conditional comparison density d(v;
e Y (t)) via
simulate from F [y; Y (t + h)|Y (t) = Q(u; Y (t))] by utilizing the given sample Q(u;
b Y (t + h)|Y (t) =
accept-reject rule to arrive at the “smooth” nonparametric model for Q[v;
Q(u; Y (t)]. See Parzen and Mukhopadhyay (2013b) for details about the method. Our
proposed algorithm generates “large” additional simulated samples from the conditional
distribution, which allows us to accurately estimate the conditional quantiles (especially
the extreme quantiles). By construction, our method guaranteed to produce non-crossing
quantile curves - thus tackles a challenging practical problem.
For S&P 500 data we first nonparametrically estimate the conditional comparison densities
b h, u) shown in Left panel of Fig 9 for F (y; Y (t)) = .01, .5 and .99, which can be thought of
d(v;
as a “weighting function” for unconditional marginal distribution to produce the conditional
distributions:




fb y; Y (t + h)|Y (t) = Q(u; Y (t)) = f (y; Y (t)) × db F (y; Y (t + h)); h, u .
14

(3.17)

f(y;Y(t+1)|Y(t)=Q(.01;X))

0.20

Density

3
2

0.00

1

0.10

4

0.30

d[v;Y(t+1),Y(t+1)|Y(t)=Q(.01;Y(t))]

0.0

0.2

0.4

0.6

0.8

−20

1.0

−15

−10

−5

0

5

10

f(y;Y(t+1)|Y(t)=Q(.5;X))

0.4

0.4
0.0

0.6

0.2

0.8

Density

1.0

0.6

1.2

d[v;Y(t+1),Y(t+1)|Y(t)=Q(.5;Y(t))]

0.0

0.2

0.4

0.6

0.8

−6

1.0

−4

−2

0

2

4

6

f(y;Y(t+1)|Y(t)=Q(.99;X))

0.3
0.2
0.0

1

0.1

2

Density

3

4

0.4

d[v;Y(t+1),Y(t+1)|Y(t)=Q(.99;Y(t))]

0.0

0.2

0.4

0.6

0.8

−5

1.0

0

5

10

Figure 9: Each row display the estimated conditional comparison density and the corresponding conditional distribution for u = .01, .5, .99.
This density estimation technique belongs to the skew-G modeling class (Mukhopadhyay,
2016). We simulate n = 10, 000 samples from the fb(y; Y (t + h)|Y (t)) by accept-reject
b h, u), u = {.01, .5, .99}. The histograms and the smooth conditional
sampling from d(v;
densities are shown in the right panel of Fig 9. It shows some typical shapes in terms of
long-tailedness.
b Y (t+h)|Y (t)), for
Next we proceed to estimate the nonparametric conditional quantiles Q(v;
v = .001, .25, .5, .75, .999 from the simulated data. Fig 10 shows the estimated conditional
quantiles. The extreme conditional quantiles have a special significance in the context
of financial time series. They sometimes popularly known as Conditional Value at Risk
(CoVaR) - currently the most popular quantitative risk management tool (see Adrian and
Brunnermeier (2011), Engle and Manganelli (2004)). The red solid line in Fig 10 is the


b .001; Y (t + 1) | Y (t) = Q(u; Y (t)) , which is known as .1% CoVar function for one-day
Q
holding period for S&P 500 daily return data. Although the upper conditional quantile curve
b
Q(.999;
Y (t + 1)|Y (t)) (blue solid line) show symmetric behaviour around F (y; Y (t)) = .5,
15

10

●
●

−5
−20

−15

−10

Y(t+1)

0

5

●
●●
●
●
●
●
●
●
●●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
● ●●
●●● ●
● ●
●
● ● ●
● ●● ●●
●
●
●
●●
●●● ● ●
●
●
●
● ● ●●●
●● ●●
● ●●
●●
●●
●● ● ●
●●● ●● ●●●
●
●
●
●
● ● ●
●●
●● ● ● ● ●●●
●●
●
●●●● ● ●
●
●
●
●● ●
●
●
●
●
●
●
● ●● ●●●●
●
●
●
●
●
● ●●
●
●
●●●
● ●●●
●
●
●
●
●
●
●●
●●
●
●
●
●
●
●●
●
●
●●
●
●●
●
●
●●
●
●
●
●
●
●
●●
●●●●●●
●●
●
●
●●
● ●●
●●
●
●
●
●●
●●
●
●
●
●
●●●●●●
●●●
● ●●●●
●
●
●
●
●
●
● ●●
●●●
●
●
●
●●
● ●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●●●
●
●
●
●●
●
●
●
●
●●
●
●●
●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●●
●
●●
●
●
●
●
●
●
●
●●●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●●
●
●
●
●●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
● ●● ●
●
●●● ●●
●
●
●●
●
●●
●
●●●●●
●
●
●
●●
●●
●
●
●●●●
●●
●
●
●●
●●●
●
●
●●
●●
●●
● ●●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●●●●●
● ●●
●
●
●●
●● ●●●
●
●
●●
●
●
●●●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●●
●
●
●●
●
●
●
●●●●
●●
●
●
●
●
●
●
●●
●
●
●
●●●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●
●
●
●
●
●●
●●
●
●
●
●
●
●
●
●
●●
●●
●●
●●
●●●●●
● ●
●
●
●
●
●
●●
●●●
●
●●●
●
●
●
●
●●●
●●●
●
●
●●● ●●
●
●●●●
●
●
●
●
●
●
●●
●●●
●
●
●
●
●
●
●●
●●
● ●●●
●
●●
●
●
●
●
● ● ●●●●●●●
●
●●
●
●●●●● ●●
●●
●
●
●●●●●●
●●●●
●● ●
●●
● ● ●●
●
●
●
●●●●●●●
● ●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●● ● ● ●●●●●
●●●
●
●
●
● ●
●●●
●● ●●
● ●
●
● ●●●
●●●
●● ● ●●●● ●● ●●
●●
●●
●●
● ●
●
●
●
●
●
●
●
●● ● ●
●●●
●
●
●●
●
●
● ●●
● ●●
●
●
●●
●●
●
●
●
●
●
●
●
● ●
●
●●
●
●
●
● ●
●
●
●
●
●
●
●
●
●
●
●
●

Q(v=.001;Y(t+1)|Y(t))
Q(v=.25;Y(t+1)|Y(t))
Q(v=.50;Y(t+1)|Y(t))
Q(v=.75;Y(t+1)|Y(t))
Q(v=.999;Y(t+1)|Y(t))

●

0.0

0.2

0.4

0.6

0.8

1.0

Fmid(y;Y(t))

Figure 10: The figure shows estimated non-parametric conditional quantile curves for S&p
b
500 return data. The red solid line which represent Q(.001;
Y (t+1)|Y (t)) is popularly known
as one day .1% conditional value at risk measure (CoVaR).
the lower quantile has a prominent asymmetric shape. These conditional quantiles give
ultimate description of the auto-regressive dependence of S&P 500 return movement in the
tail-region.

3.6

Nonlinear Spectrum Analysis

Here we extends the concept of spectral density for nonlinear processes. Display the LPSpectrum - autoregressive (AR) spectral density estimates of YS1 (t), . . . , YS4 (t). Spectral
density for each LP-transformed series is defined as
X
f (ω; YSj ) =
LP[j, j; Y (t), Y (t + h)] e−i2πhω ,

−1/2 < ω < 1/2

h

=

X

Cov[YSj (t), YSj (t + h)] e−i2πhω ,

−1/2 < ω < 1/2.

(3.18)

h

We separately fit the univariate AR model for the components of Vec(YS)(t) and use BIC order selection criterion to select the “best” parsimonious parametrization using Burg method.
Finally we use the estimated model coefficients to produce the “smooth” estimate of spectral
density function (see Eq. 1.2). The copula spectral density is defined as
X
f (ω; u, v) =
cop(u, v; h) e−i2πhω ,
−1/2 < ω < 1/2.

(3.19)

h

To estimate the copula spectral density we use the LP comoment based nonparametric
copula density estimate. Note that both the serial copula (3.12) and the corresponding
16

spectral density (3.25) captures the same amount of information for serial dependence of
{Y (t)}. For that reason, we recommend to compute AutoLPinfor as a general dependence
measure for non-gaussian nonlinear processes. Application of our LPSpectral tool on S&P
LPSpectrum

2
1
0

Log Spectra

3

Raw data (order=2)
S1 (order=1)
S2 (order=25)
S3 (order=0)
S4 (order=13)

0.0

0.1

0.2

0.3

0.4

0.5

Frequency

Figure 11: LPSpec: AR spectral density estimate for S&P 500 return data. Order selected
by BIC method. This provides a diagnostic tool for providing evidence of hidden periodicities
in non-Gaussian nonlinear time series.
500 return data is shown in Fig 11. Few interesting observations: (i) the conventional
spectral density (black solid line) provides no insight into the (complex) serial dependency
present in the data. (ii) The nonlinearity in the series is captured by the interesting shapes
of our specially-designed times series YS2 (t) and YS4 (t), which classical (linear) correlogram
based spectra can not account for. (iii) The shape of the spectra of Z(Y (t)) and the ranktransformed time series YS1 (t) looks very similar. and a (iv) pronounced singularity near
zero of the spectrum of YS2 (t) hints some kind of “long memory” behavior. This phenomena
also known as regular variation representation at frequency ω = 0 (Granger and Joyeux,
1980).
A quick diagnostic measure for screening significant spectrums can be computed via the
R 1/2
information number 2 0 log fb(ω; Sj ) dω. The LPSpectrum methodology is highly robust
and thus can tackle the heavy-tailed S&P data quite successfully.

17

3.7

Nonparametric Model Specification

The ultimate goal of empirical time series analysis is nonparametric model identification. To
model the univariate stationary nonlinear process we specify multiple autoregressive model

T
based on Vec(YS)(t) = YS1 (t), . . . , YSk (t) of the form
Vec(YS)(t) =

m
X

A(k; m) Vec(YS)(t − k) + (t).

(3.20)

k=1

where (t) is multivariate mean zero Gaussian white noise with covariance Σm . These system
of equation jointly describe the dynamics of the nonlinear process and how it evolves over
time. We use BIC criterion to select the model order which selects the model order m which
minimizes

b m | + mk2 log T .
BIC(m) = log |Σ
T

(3.21)

We carry out this steps for our S&P 500 return data. We estimate our multiple AR model

T
based on Vec(YS)(t) = YS1 (t), YS2 (t), YS4 (t) . We discard YS3 (t) due to it’s flat spectrum (see Fig 11). BIC selects “best” order 8. Although, the complete description of the
estimated model is clearly cumbersome, we provide below the approximate structure by
selecting few large coefficients from the actual matrix equation. The goal is to interpret
the coefficients (statistical parameters) of the estimated model and relate it with economic
theory (scientific parameters/theory). This multiple AR LP-model (in short we call LPVAR)
is given by
YS1 (t) ≈ .071 YS1 (t − 1) − .024 YS1 (t − 2) + 1 (t)
YS2 (t) ≈ −.063 YS1 (t − 1) − .075 YS1 (t − 2) + .06 YS2 (t − 2) + .123 YS2 (t − 5) + .04 YS4 (t − 2) + 2 (t)
YS4 (t) ≈ .04 YS4 (t − 1) + .038 YS4 (t − 2) + .04 YS2 (t − 3) + 4 (t).

and the residual covariance matrix is


0.993 −0.001 −0.002

c8 = 
Σ
−0.001 0.853 −0.058
−0.002 −0.058 0.964
The autoregressive model of YS2 (t) can be considered as a robust stock return volatility
model - LPVolatility Modeling, which is less affected by unusually large extreme events.
The model for YS2 (t) automatically discover many known facts: (a) the sign of the coefficient linking volatility and return is negative -confirming “leverage effect”; (b) YS2 (t) is
positively autocorrelated, known as volatility clustering; (c) The positive interaction with
lagged YS4 (t) accounts for the “excess kurtosis”.

4

Conclusion

This article provides a pragmatic and comprehensive framework for nonlinear time series
modeling that is easier to use, more versatile and has a strong theoretical foundation based
18

(3.22)

on recently developed theory on unified algorithms of data science via LP modeling (Parzen
and Mukhopadhyay, 2012, 2013a,b, Mukhopadhyay and Parzen, 2014, Mukhopadhyay, 2016,
2017). Summary and broader implications of the proposed research:
• From the theoretical standpoint, the unique aspect of our proposal lies in its ability to
simultaneously embrace and employ spectral domain, time domain, quantile domain, and
the information domain analysis for enhanced insights, which to the best of our knowledge
has not appeared in the nonlinear time series literature before.
• From a practical angle, the novelty of our technique is that it permits us to use the
techniques from linear Gaussian time series to create non-Gaussian nonlinear time series
models with highly interpretable parameters. This aspect makes LPTime computationally
extremely attractive for the data scientists as they can now borrow all the standard time
series analysis machinery from R libraries for implementation purpose.
• From the pedagogical side, we believe that these concepts and methods can easily be
augmented with the standard time series analysis course to modernize the current curriculum
so that students can handle complex time series modeling problems (McNeil et al., 2010)
using the tools they are already familiar with.
The main thrust of this article is to describe and interpret the steps of LPTime technology
to create a realistic general-purpose algorithm for empirical time series modeling. In addition, many new theoretical results and diagnostic measures are presented which laid the
foundation for the algorithmic implementation of LPTime. We showed how LPTime can
systematically explore the data to discover empirical facts hidden in time series. For example, LPTime empirical modeling of S&P 500 return data reproduces the ‘stylized facts’:
(a) heavy tails; (b) non-Gaussian;(c) nonlinear serial dependence; (d) tail-correlation; (e)
asymmetric dependence; (f) volatility clustering; (g) long-memory volatility structure; (h)
efficient market hypothesis; (i) leverage effect; (j) excess kurtosis, in a coherent manner
under single general unified framework. We have emphasized how the statistical parameters
of our model can be interpreted in the light of established economic theory.
We have recently applied this theory for large-scale eye-movement pattern discovery problem,
which came out as the winner (among 82 competing algorithms) of the 2014 IEEE International Biometric Eye Movements Verification and Identification Competition (Mukhopadhyay and Nandi, 2017).
We conclude with some general references. Few popular articles Granger (2003, 1993),
Granger and Lin (1994), Engle (1982), Parzen (1979, 1967), Tukey (1980), Brillinger (1977,
2004), Salmon (2012); Books Terasvirta et al. (2010), Woodward et al. (2011), Tsay (2010);
and Review articles Granger (1998), Hendry (2011).
The proposed algorithm is implemented in the R package LPTime, which is available on CRAN.
19

References
Adrian, T. and M. K. Brunnermeier (2011). Covar. Technical report, National Bureau of
Economic Research.
Brillinger, D. R. (1977). The identification of a particular nonlinear time series system.
Biometrika, 509–515.
Brillinger, D. R. (2004). Some data analyses using mutual information. Brazilian Journal
of Probability and Statistics 18 (6), 163–183.
Engle, R. F. (1982). Autoregressive conditional heteroscedasticity with estimates of the
variance of united kingdom inflation. Econometrica: Journal of the Econometric Society,
987–1007.
Engle, R. F. and S. Manganelli (2004). Caviar: Conditional autoregressive value at risk by
regression quantiles. Journal of Business & Economic Statistics 22 (4), 367–381.
Granger, C. and J.-L. Lin (1994). Using the mutual information coefficient to identify lags
in nonlinear models. Journal of time series analysis 15 (4), 371–384.
Granger, C. W. (1993). Strategies for modelling nonlinear time-series relationships. Economic Record 69 (3), 233–238.
Granger, C. W. (1998). Overview of nonlinear time series specification in economics. NSF
Symposium on Nonlinear Time Series Models, University of California, Berkeley.
Granger, C. W. (2003). Time series concepts for conditional distributions. Oxford Bulletin
of Economics and Statistics 65 (s1), 689–701.
Granger, C. W. and R. Joyeux (1980). An introduction to long-memory time series models
and fractional differencing. Journal of time series analysis 1 (1), 15–29.
Hendry, D. F. (2011). Empirical economic model discovery and theory evaluation. Rationality, Markets and Morals 2, 115–145.
McNeil, A. J., R. Frey, and P. Embrechts (2010). Quantitative risk management: concepts,
techniques, and tools. Princeton university press.
Mukhopadhyay, S. (2016). Large scale signal detection: A unifying view. Biometrics 72 (2),
325–334.
Mukhopadhyay, S. (2017). Large-scale mode identification and data-driven sciences. Electronic Journal of Statistics 11 (1), 215–240.
20

Mukhopadhyay, S. and S. Nandi (2017). LPiTrack: Eye movement pattern recognition algorithm and application to biometric identification. Machine Learning, doi 10.1007/s10994017-5649-1 .
Mukhopadhyay, S. and E. Parzen (2014). LP approach to statistical modeling. Preprint
arXiv:1405.2601 .
Parzen, E. (1967). On empirical multiple time series analysis. In Proc. Fifth Berkeley Symp.
Math. Statist. Prob, Volume 1, pp. 305–340.
Parzen, E. (1979). Nonparametric statistical data modeling (with discussion). Journal of
the American Statistical Association 74, 105–131.
Parzen, E. (1992). Time series, statistics, and information. New Directions in Time Series
Analysis (ed. E. Parzen et al). Springer Verlag: New York., 265–286.
Parzen, E. (1997). Comparison distributions and quantile limit theorems. International
Conference on Asymptotic Methods in Probability and Statistics (Carleton University,
Canada, 8-13 July 1997).
Parzen, E. and S. Mukhopadhyay (2012). Modeling, Dependence, Classification, United
Statistical Science, Many Cultures. arXiv:1204.4699 .
Parzen, E. and S. Mukhopadhyay (2013a). United Statistical Algorithms, LP comoment,
Copula Density, Nonparametric Modeling. 59th ISI World Statistics Congress (WSC),
Hong Kong.
Parzen, E. and S. Mukhopadhyay (2013b). United Statistical Algorithms, Small and Big
Data, Future of Statisticians. arXiv:1308.0641 .
Salmon, F. (2012). The formula that killed wall street. Significance 9 (1), 16–20.
Terasvirta, T., D. Tjostheim, and C. W. Granger (2010). Modelling nonlinear economic
time series. OUP Catalogue.
Tsay, R. S. (2010). Analysis of financial time series. Wiley.
Tukey, J. (1980). Can we predict where “time series” should go next. Directions in Times
Series, IMS, Hayward, Calif , 1–31.
Woodward, W. A., H. L. Gray, and A. C. Elliott (2011). Applied time series analysis. CRC
press.

21

