Merge aa75ba3 into 173ecde

.gitignore

@@ -5,7 +5,7 @@
 *.log
 docs/build/
 docs/site
-src/data
+src/data/bollerslev_ghysels.txt
 docs/Manifest.toml
 Manifest.toml
 make-require.jl

.travis.yml

@@ -35,5 +35,5 @@ jobs:
                                     Pkg.develop(PackageSpec(path=pwd()));
                                     Pkg.instantiate();
                                     Pkg.build("ARCHModels")' #https://github.com/JuliaLang/Pkg.jl/issues/862
-        - julia --project=docs/ docs/make.jl
+        - julia --depwarn=no --project=docs/ docs/make.jl
       after_success: skip

deps/build.jl

@@ -1,5 +1,28 @@
+#using Retry
+using DelimitedFiles
 datadir = joinpath(@__DIR__, "..", "src", "data")
-if ~isdir(datadir)
-    mkdir(datadir)
-end
-download("http://people.stern.nyu.edu/wgreene/Text/Edition7/TableF20-1.txt", joinpath(datadir, "bollerslev_ghysels.txt"))
+isdir(datadir) || mkdir(datadir)
+@info "Downloading Bollerslev and Ghysels data..."
+isfile(joinpath(datadir, "bollerslev_ghysels.txt")) || download("http://people.stern.nyu.edu/wgreene/Text/Edition7/TableF20-1.txt", joinpath(datadir, "bollerslev_ghysels.txt"))
+
+# @info "Downloading stock data..."
+# #"DOW" is excluded because it's listed too late
+# tickers = ["AAPL", "IBM", "XOM", "KO", "MSFT", "INTC", "MRK", "PG", "VZ", "WBA", "V", "JNJ", "PFE", "CSCO", "TRV", "WMT", "MMM", "UTX", "UNH", "NKE", "HD", "BA", "AXP", "MCD", "CAT", "GS", "JPM", "CVX", "DIS"]
+# alldata = zeros(2786, 29)
+# for (j, ticker) in enumerate(tickers)
+#         @repeat 4 try
+#             @info "...$ticker"
+#             filename = joinpath(datadir, "$ticker.csv")
+#             isfile(joinpath(datadir, "$ticker.csv")) || download("http://quotes.wsj.com/$ticker/historical-prices/download?num_rows=100000000&range_days=100000000&startDate=03/19/2008&endDate=04/11/2019", filename)
+#             data = parse.(Float64, readdlm(joinpath(datadir, "$ticker.csv"), ',', String, skipstart=1)[:, 5])
+#             length(data) == 2786 || error("Download failed for $ticker.")
+#             alldata[:, j] .= data
+#             rm(filename)
+#         catch e
+#             @delay_retry if 1==1 end
+#         end
+# end
+# alldata = 100 * diff(log.(alldata), dims=1)
+# open(joinpath(datadir, "dow29.csv"), "w") do io
+#      writedlm(io, alldata, ',')
+# end

docs/make.jl

@@ -7,8 +7,8 @@ makedocs(modules=[ARCHModels],
         doctest=true,
         strict=true,
         pages = ["Home" => "index.md",
-                 "types.md",
-                 "manual.md",
+                 "Univariate ARCH Models" => Any["univariateintro.md", "univariateusage.md"],
+                 "Multivariate ARCH Models" => Any["multivariateintro.md", "multivariateusage.md"],
                  "reference.md"
                  ]
         )

docs/src/index.md

@@ -33,12 +33,6 @@ add ARCHModels
 
 in the Pkg REPL mode (which is entered by pressing `]` at the prompt).
 
-## Contents
-```@contents
-Pages = ["types.md", "manual.md", "reference.md"]
-Depth = 2
-```
-
 ## Acknowledgements
 
 This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement No. 750559).

docs/src/multivariateintro.md

@@ -0,0 +1,41 @@
+# Introduction
+
+ARCH models naturally generalize to the multivariate setting. Consider a time series of daily asset returns $\{r_t\}_{t\in 1, \ldots, T}$, where now  $r_t\in\mathbb{R}^d$. Similarly to the univariate case, the general model structure is
+
+\[
+r_t=\mu_t+\Sigma_t^{1/2}z_t,\quad \mu_t\equiv\mathbb{E}[r_t\mid\mathcal{F}_{t-1}],\quad \Sigma_t\equiv\mathbb{E}[(r_t-\mu_t)(r_t-\mu_t)^\mathrm{T}]\mid\mathcal{F}_{t-1}].
+\]
+
+A multivariate ARCH model specifies the conditional covariance matrix $\Sigma_t$ in terms of past returns, conditional (co)variances, and potentially other variables. The main challenge in multivariate ARCH modelling is the _curse of dimensionality_: allowing each of the $(d)(d+1)/2$ elements of $\Sigma_t$ to depend on the past returns of all $d$ other assets requires $O(d^4)$ parameters without imposing additional structure. Multivariate ARCH models differ in which structure they impose.
+
+The following multivariate ARCH models are currently available in this package:
+
+  * The CCC model of [Bollerslev (1990)](https://doi.org/10.2307/2109358)
+  * The DCC model of [Engle (2002)](https://doi.org/10.1198/073500102288618487)
+
+These may be combined with different mean specifications as in the univariate case, and (in principle) with different specifications for the joint distribution of the $z_t$, although at present, only the multivariate standard normal is supported. Multivariate ARCH models are represented as instances of [`MultivariateARCHModel`](@ref), which like [`UnivariateARCHModel`](@ref) subtypes [`ARCHModel`](@ref).
+
+# Type hierarchy
+## [Covariance specifications](@id covspec)
+Volatility specifications describe the evolution of $\Sigma_t$. They are modelled as subtypes of [`MultivariateVolatilitySpec`](@ref).
+
+### CCC
+The CCC (and DCC, see below) models are examples of _conditional correlation_ models. They decompose
+$\Sigma_t$ as
+\[
+\Sigma_t=D_t R_t D_t,
+\]
+where $R_t$ is the conditional correlation matrix and $D_t$ is a diagonal matrix containing the volatilities of the individual assets, which are modelled as univariate ARCH processes. In the constant conditional correlation (CCC) model, $R_t=R$ is assumed constant. The model is typically, including here, estimated in a two-step procedure: first, univariate ARCH models are fitted to the $d$ asset returns, and then $R$ is estimated as the sample correlation matrix of the standardized residuals.
+
+### DCC
+The DCC model extends the CCC model by making the $R_t$ dynamic (hence the name, dynamic conditional correlation model). In particular,
+
+\[
+R_{ij, t} = \frac{Q_{ij,t}}{\sqrt{Q_{ii,t}Q_{jj,t}}}, \
+\]
+where
+\[Q_{t} =(1-\theta _{1}-\theta _{2})\bar{Q}+\theta _{1}\epsilon
+_{t-1}\epsilon _{t-1}^{\prime }+\theta _{2}Q_{t-1},
+\]
+$\epsilon _{t}=D_{t}^{-1}a_{t}$, $Q_{t}=\mathrm{cov}%
+(\epsilon _{t}|F_{t-1})$, and $\bar{Q}=\mathrm{cov}(\epsilon _{t})$.

docs/src/multivariateusage.md

@@ -0,0 +1 @@
+# Usage

docs/src/types.md → docs/src/univariateintro.md

@@ -1,4 +1,4 @@
-# Introduction and type hierarchy
+# Introduction
 Consider a sample of daily asset returns ``\{r_t\}_{t\in\{1,\ldots,T\}}``. All models covered in this package share the same basic structure, in that they decompose the return into a conditional mean and a mean-zero innovation:
 ```math
 r_t=\mu_t+\sigma_tz_t,\quad \mu_t\equiv\mathbb{E}[r_t\mid\mathcal{F}_{t-1}],\quad \sigma_t^2\equiv\mathbb{E}[(r_t-\mu_t)^2\mid\mathcal{F}_{t-1}],
@@ -7,8 +7,9 @@ where ``z_t`` is identically and independently distributed according to some law
 
 This package represents a univariate (G)ARCH model as an instance of [`UnivariateARCHModel`](@ref), which implements the interface of `StatisticalModel` from [`StatsBase`](http://juliastats.github.io/StatsBase.jl/stable/statmodels.html). An instance of this type contains a vector of data (such as equity returns), and encapsulates information about the [volatility specification](@ref volaspec) (e.g., [GARCH](@ref) or [EGARCH](@ref)), the [mean specification](@ref meanspec) (e.g., whether an intercept is included), and the [error distribution](@ref Distributions).
 
+# Type hierarchy
 ## [Volatility specifications](@id volaspec)
-Volatility specifications describe the evolution of ``\sigma_t``. They are modelled as subtypes of [`VolatilitySpec`](@ref). There is one type for each class of (G)ARCH model, parameterized by the number(s) of lags (e.g., ``p``, ``q`` for a GARCH(p, q) model). For each volatility specification, the order of the parameters in the coefficient vector is such that all parameters pertaining to the first type parameter (``p``) appear before those pertaining to the second (``q``).
+Volatility specifications describe the evolution of ``\sigma_t``. They are modelled as subtypes of [`UnivariateVolatilitySpec`](@ref). There is one type for each class of (G)ARCH model, parameterized by the number(s) of lags (e.g., ``p``, ``q`` for a GARCH(p, q) model). For each volatility specification, the order of the parameters in the coefficient vector is such that all parameters pertaining to the first type parameter (``p``) appear before those pertaining to the second (``q``).
 ### ARCH
 With ``a_t\equiv r_t-\mu_t``, the ARCH(q) volatility specification, due to [Engle (1982)](https://doi.org/10.2307/1912773 ), is
 ```math
@@ -160,7 +161,7 @@ julia> const MyStdT = Standardized{TDist};
 julia> ARCHModels.startingvals(::Type{<:MyStdT}, data::Vector{T}) where T = T[3.]
 ```
 ## Working with UnivariateARCHModels
-The constructor for [`UnivariateARCHModel`](@ref) takes two mandatory arguments: an instance of a subtype of [`VolatilitySpec`](@ref), and a vector of returns. The mean specification and error distribution can be changed via the keyword arguments `meanspec` and `dist`, which respectively default to `NoIntercept` and `StdNormal`.
+The constructor for [`UnivariateARCHModel`](@ref) takes two mandatory arguments: an instance of a subtype of [`UnivariateVolatilitySpec`](@ref), and a vector of returns. The mean specification and error distribution can be changed via the keyword arguments `meanspec` and `dist`, which respectively default to `NoIntercept` and `StdNormal`.
 
 For example, to construct a GARCH(1, 1) model with an intercept and ``t``-distributed errors, one would do
 ```jldoctest TYPES
@@ -232,3 +233,11 @@ julia> nobs(am)
 ```
 
 Other useful methods include [`means`](@ref), [`volatilities`](@ref) and [`residuals`](@ref).
+
+```@meta
+DocTestSetup = quote
+    using Random
+    Random.seed!(1)
+end
+DocTestFilters = r".*[0-9\.]"
+```

docs/src/manual.md → docs/src/univariateusage.md

@@ -1,10 +1,3 @@
-```@meta
-DocTestSetup = quote
-    using Random
-    Random.seed!(1)
-end
-DocTestFilters = r".*[0-9\.]"
-```
 # Usage
 We will be using the data from [Bollerslev and Ghysels (1986)](https://doi.org/10.2307/1392425), available as the constant [`BG96`](@ref). The data consist of daily German mark/British pound exchange rates (1974 observations) and are often used in evaluating
 implementations of (G)ARCH models (see, e.g., [Brooks et.al. (2001)](https://doi.org/10.1016/S0169-2070(00)00070-4). We begin by convincing ourselves that the data exhibit ARCH effects; a quick and dirty way of doing this is to look at the sample autocorrelation function of the squared returns:
@@ -81,7 +74,7 @@ This returns an instance of [`UnivariateARCHModel`](@ref), as described in the s
 
 The [`fit`](@ref) method supports a number of keyword arguments; the full signature is
 ```julia
-fit(::Type{<:VolatilitySpec}, data::Vector; dist=StdNormal, meanspec=Intercept, algorithm=BFGS(), autodiff=:forward, kwargs...)
+fit(::Type{<:UnivariateVolatilitySpec}, data::Vector; dist=StdNormal, meanspec=Intercept, algorithm=BFGS(), autodiff=:forward, kwargs...)
 ```
 
 Their meaning is as follows:
@@ -141,7 +134,7 @@ Distribution parameters:
 ─────────────────────────────────────────
 ```
 
-An alternative approach to fitting a [`VolatilitySpec`](@ref) to `BG96` is to first construct
+An alternative approach to fitting a [`UnivariateVolatilitySpec`](@ref) to `BG96` is to first construct
 a [`UnivariateARCHModel`](@ref) containing the data, and then using [`fit!`](@ref) to modify it in place:
 
 ```jldoctest MANUAL
@@ -231,7 +224,7 @@ Volatility parameters:
 ```
 ## Model selection
 The function [`selectmodel`](@ref) can be used for automatic model selection, based on an information crititerion. Given
-a class of model (i.e., a subtype of [`VolatilitySpec`](@ref)), it will return a fitted [`UnivariateARCHModel`](@ref), with the lag length
+a class of model (i.e., a subtype of [`UnivariateVolatilitySpec`](@ref)), it will return a fitted [`UnivariateARCHModel`](@ref), with the lag length
 parameters (i.e., ``p`` and ``q`` in the case of [`GARCH`](@ref)) chosen to minimize the desired criterion. The [BIC](https://en.wikipedia.org/wiki/Bayesian_information_criterion) is used by default.
 
 Eg., the following selects the optimal (minimum AIC) EGARCH(o, p, q) model, where o, p, q < 2,  assuming ``t`` distributed errors.
@@ -410,7 +403,7 @@ Details:
 ```
 By default, the number of lags is chosen as the maximum order of the volatility specification (e.g., ``\max(p, q)`` for a GARCH(p, q) model). Here, the test does not reject, indicating that a GARCH(1, 1) specification is sufficient for modelling the volatility clustering (a common finding).
 ## Simulation
-To simulate from a [`UnivariateARCHModel`](@ref), use [`simulate`](@ref). You can either specify the [`VolatilitySpec`](@ref) (and optionally the distribution and mean specification) and desired number of observations, or pass an existing [`UnivariateARCHModel`](@ref). Use [`simulate!`](@ref) to modify the data in place. Example:
+To simulate from a [`UnivariateARCHModel`](@ref), use [`simulate`](@ref). You can either specify the [`UnivariateVolatilitySpec`](@ref) (and optionally the distribution and mean specification) and desired number of observations, or pass an existing [`UnivariateARCHModel`](@ref). Use [`simulate!`](@ref) to modify the data in place. Example:
 ```jldoctest MANUAL
 julia> am3 = simulate(GARCH{1, 1}([1., .9, .05]), 1000; warmup=500, meanspec=Intercept(5.), dist=StdT(3.))
 

src/ARCHModels.jl

@@ -9,7 +9,7 @@
 #a simulated AM should probably contain a (zero) intercept, so that fit! is consistent with fit.
 #the constructor for UnivariateARCHModel should make a copy of its args
 #implement lrtest
-#allow uninititalized constructors for VolatilitySpec, MeanSpec and StandardizedDistribution? If so, then be consistent with how they are defined
+#allow uninititalized constructors for UnivariateVolatilitySpec, MeanSpec and StandardizedDistribution? If so, then be consistent with how they are defined
 #  (change for meanspec and dist ), document, and test. Also, NaN is prob. safer than undef.
 #logconst needs to return the correct type
 """
@@ -29,6 +29,7 @@ using Roots
 using LinearAlgebra
 using DataStructures: CircularBuffer
 using DelimitedFiles
+using Statistics: cov
 
 import Distributions: quantile
 import Base: show, showerror, eltype
@@ -38,10 +39,12 @@ import HypothesisTests: HypothesisTest, testname, population_param_of_interest,
 import StatsBase: StatisticalModel, stderror, loglikelihood, nobs, fit, fit!, confint, aic,
                   bic, aicc, dof, coef, coefnames, coeftable, CoefTable,
 				  informationmatrix, islinear, score, vcov, residuals, predict
-export ARCHModel, UnivariateARCHModel, VolatilitySpec, StandardizedDistribution, Standardized, MeanSpec,
+export ARCHModel, UnivariateARCHModel, UnivariateVolatilitySpec, StandardizedDistribution, Standardized, MeanSpec,
        simulate, simulate!, selectmodel, StdNormal, StdT, StdGED, Intercept, Regression,
-       NoIntercept, ARMA, AR, MA, BG96, volatilities, mean, quantile, VaRs, pvalue, means,
-	   EGARCH, ARCH, GARCH, TGARCH, ARCHLMTest, DQTest
+       NoIntercept, ARMA, AR, MA, BG96, volatilities, mean, quantile, VaRs, pvalue, means, VolatilitySpec,
+	   MultivariateVolatilitySpec, MultivariateStandardizedDistribution, MultivariateARCHModel, MultivariateStdNormal,
+	   EGARCH, ARCH, GARCH, TGARCH, ARCHLMTest, DQTest,
+	   DOW29, DCC, CCC, covariances, correlations
 
 
 include("utils.jl")
@@ -52,14 +55,17 @@ include("univariatestandardizeddistributions.jl")
 include("EGARCH.jl")
 include("TGARCH.jl")
 include("tests.jl")
+include("multivariatearchmodel.jl")
+include("multivariatestandardizeddistributions.jl")
+include("DCC.jl")
 function __init__()
 	@require GLM = "38e38edf-8417-5370-95a0-9cbb8c7f171a" begin
 		using .GLM
 		import .StatsModels: TableRegressionModel
-		function fit(vs::Type{VS}, lm::TableRegressionModel{<:LinearModel}; kwargs...) where VS<:VolatilitySpec
+		function fit(vs::Type{VS}, lm::TableRegressionModel{<:LinearModel}; kwargs...) where VS<:UnivariateVolatilitySpec
 			fit(vs, response(lm.model); meanspec=Regression(modelmatrix(lm.model); coefnames=coefnames(lm)), kwargs...)
 		end
-		function fit(vs::Type{VS}, lm::LinearModel; kwargs...) where VS<:VolatilitySpec
+		function fit(vs::Type{VS}, lm::LinearModel; kwargs...) where VS<:UnivariateVolatilitySpec
 			fit(vs, response(lm); meanspec=Regression(modelmatrix(lm)), kwargs...)
 		end
 	end