Initial Commit

Author: CHNPAT005

.DS_Store

@@ -0,0 +1,93 @@
+## Author: Patrick Chang
+# Script file for the MM Complex Fourier Transform (Dirichlet)
+# Supporting Algorithms are at the start of the script
+#  Include:
+#           - Scale function to re-scale time to [0, 2 \pi]
+# Number of Fourier Coefficients automatically chosen so that events
+# are not aliased
+
+#---------------------------------------------------------------------------
+
+### Data Format:
+## p = [n x m] matrix of prices, log returns are computed in the function
+# non-trading times are indicated by NaNs
+## t = [n x m] matrix of trading times, non-trading times are indicated by NaNs
+# dimensions of p and t must match.
+## N = Optional input for cutoff frequency
+
+#---------------------------------------------------------------------------
+
+using ArgCheck; using LinearAlgebra
+
+#---------------------------------------------------------------------------
+### Supporting functions
+
+function scale(t)
+    maxt = maximum(filter(!isnan, t))
+    mint = minimum(filter(!isnan, t))
+
+    tau = (2*pi) .* (t .- mint) ./ (maxt - mint)
+    return tau
+end
+
+#---------------------------------------------------------------------------
+# Dirichlet Kernel implementaion
+# wave range from -N:N
+# Exploits c(-k) = \bar{c(k)} in the computation so that
+# we need only compute 1:N and sum(DiffP) for efficiency
+
+function CFTcorrDK(p, t; kwargs...)
+    ## Pre-allocate arrays and check Data
+    np = size(p)[1]
+    mp = size(p)[2]
+    nt = size(t)[1]
+
+    @argcheck size(p) == size(t) DimensionMismatch
+
+    # Re-scale trading times
+    tau = scale(t)
+    # Computing minimum time change
+    # minumum step size to avoid smoothing
+    dtau = diff(filter(!isnan, tau))
+    taumin = minimum(filter((x) -> x>0, dtau))
+    taumax = 2*pi
+    # Sampling Freq.
+    N0 = taumax/taumin
+
+    # Optional Cutoff - if not specified we use Nyquist Cutoff
+    kwargs = Dict(kwargs)
+
+    if haskey(kwargs, :N)
+        k = collect(1:1:kwargs[:N])
+    else
+        k = collect(1:1:floor(N0/2))
+    end
+
+    Den = length(k)
+
+    #------------------------------------------------------
+
+    c_pos = zeros(ComplexF64, mp, 2*Den + 1)
+    c_neg = zeros(ComplexF64, mp, 2*Den + 1)
+
+    for i in 1:mp
+        psii = findall(!isnan, p[:,i])
+        P = p[psii, i]
+        Time = tau[psii, i]
+        DiffP = diff(log.(P))
+
+        C = DiffP' * exp.(-1im * Time[1:(end-1),:] * k')
+
+        c_pos[i,:] = [C sum(DiffP) conj(C)]
+        c_neg[i,:] = [conj(C) sum(DiffP) C]
+    end
+
+    Sigma = 0.5 / (2*Den + 1) .* (c_pos*c_pos' + c_neg*c_neg')
+
+    Sigma = real(Sigma)
+    var = diag(Sigma)
+    sigma = sqrt.(var)
+    rho = Sigma ./ (sigma * sigma')
+
+    return rho, Sigma, var
+end

Computed Data/.DS_Store

@@ -0,0 +1,104 @@
+## Author: Patrick Chang
+# Script file for the MM Fast Fourier Transform with Zero Padding
+# Supporting Algorithms are at the start of the script
+#  Include:
+#           - Scale function to re-scale time to [0, 2 \pi]
+#           - Zero-padding function
+# Number of Fourier Coefficients used is length of data
+
+#---------------------------------------------------------------------------
+
+### Data Format:
+## p = [n x m] matrix of prices, log returns are computed in the function
+# non-trading times are indicated by NaNs
+## t = [n x m] matrix of trading times, non-trading times are indicated by NaNs
+# dimensions of p and t must match.
+
+#---------------------------------------------------------------------------
+
+using ArgCheck; using LinearAlgebra;
+using FFTW
+
+#---------------------------------------------------------------------------
+### Supporting functions
+
+function scale(t)
+    maxt = maximum(filter(!isnan, t))
+    mint = minimum(filter(!isnan, t))
+
+    tau = (2*pi) .* (t .- mint) ./ (maxt - mint)
+    return tau
+end
+
+function zero_pad(cj, xj, N0)
+    res = zeros(N0, 1)
+    nj = length(xj)
+    for j in 1:nj
+        pos = Int(round(xj[j] * (N0/(2*pi))))
+        res[pos+1] = cj[j]
+    end
+    return res
+end
+#ZPP = zero_pad(DiffP, tau[psii, 1], Int(ceil(N0)))
+#---------------------------------------------------------------------------
+
+# Fast Fourier Transform (with Zero Padding) implementaion of the Dirichlet Kernel
+# Can be used for asynchronous data
+# No option for optional cutoff (N)
+
+function FFTZPcorrDK(p, t)
+    ## Pre-allocate arrays and check Data
+    np = size(p)[1]
+    mp = size(p)[2]
+    nt = size(t)[1]
+
+    #@argcheck size(p) == size(t) DimensionMismatch
+
+    # Re-scale trading times
+    tau = scale(t)
+    # Computing minimum time change
+    # minumum step size to avoid smoothing
+    dtau = diff(filter(!isnan, tau))
+    taumin = minimum(filter((x) -> x>0, dtau))
+    taumax = 2*pi
+    # Sampling Freq.
+    N0 = taumax/taumin
+
+    k = collect(-floor(N0/2):1:floor(N0/2))
+
+    Den = length(k)
+
+    #------------------------------------------------------
+    c_pos = zeros(ComplexF64, mp, Den)
+    c_neg = zeros(ComplexF64, mp, Den)
+
+    for i in 1:mp
+        psii = findall(!isnan, p[:,i])
+        P = p[psii, i]
+        DiffP = diff(log.(P))
+        psii = psii[1:(end-1)]
+        Time = tau[psii, i]
+        ZP = zeros(Den, 1)
+
+        ZP[1:Int(floor(N0))] = zero_pad(DiffP, Time, Int(floor(N0)))
+
+        C = fft(ZP)
+
+        C_pos = C
+        C_neg = conj(C)
+
+        c_pos[i,:] = C_pos
+        c_neg[i,:] = C_neg
+    end
+
+    Sigma = zeros(ComplexF64, mp, mp)
+
+    Sigma = 0.5 / Den .* (c_pos*c_pos' + c_neg*c_neg')
+
+    Sigma = real(Sigma)
+    var = diag(Sigma)
+    sigma = sqrt.(var)
+    rho = Sigma ./ (sigma * sigma')
+
+    return rho, Sigma, var
+end

Computed Data/EppsCorrection/GBMwExpSaples.jld

@@ -0,0 +1,77 @@
+## Author: Patrick Chang
+# Script file for the MM Fast Fourier Transform
+# Supporting Algorithms are at the start of the script
+#  Include:
+#           - Scale function to re-scale time to [0, 2 \pi]
+# Number of Fourier Coefficients used is length of data
+
+#---------------------------------------------------------------------------
+
+### Data Format:
+## p = [n x m] matrix of prices, log returns are computed in the function
+# non-trading times are indicated by NaNs
+
+#---------------------------------------------------------------------------
+
+using ArgCheck; using LinearAlgebra;
+using FFTW
+
+#---------------------------------------------------------------------------
+
+# Fast Fourier Transform implementaion of the Dirichlet Kernel
+# Can ONLY be used for fully synchronous data
+# No option for optional cutoff (N)
+
+function FFTcorrDK(p)
+    ## Pre-allocate arrays and check Data
+    np = size(p)[1]
+    mp = size(p)[2]
+    #nt = size(t)[1]
+
+    #@argcheck size(p) == size(t) DimensionMismatch
+
+    # # Re-scale trading times
+    # tau = scale(t)
+    # # Computing minimum time change
+    # # minumum step size to avoid smoothing
+    # dtau = diff(filter(!isnan, tau))
+    # taumin = minimum(filter((x) -> x>0, dtau))
+    # taumax = 2*pi
+    # # Sampling Freq.
+    # N0 = taumax/taumin
+
+    if (iseven(np))
+        k = collect(0:1:(np-1)/2)
+    else
+        k = collect(0:1:(np+1)/2)
+    end
+
+    Den = length(k)
+
+    #------------------------------------------------------
+    c_pos = zeros(ComplexF64, mp, Den + (Den-1))
+    c_neg = zeros(ComplexF64, mp, Den + (Den-1))
+
+    for i in 1:mp
+        DiffP = diff(log.(p[:,i]))
+
+        C = rfft(DiffP)
+
+        C_pos = C
+        C_neg = conj(C)
+
+        c_pos[i,:] = [C_pos; C_neg[2:end]]
+        c_neg[i,:] = [C_neg; C_pos[2:end]]
+    end
+
+    Sigma = zeros(ComplexF64, mp, mp)
+
+    Sigma = 0.5 / (Den + (Den-1)) .* (c_pos*c_pos' + c_neg*c_neg')
+
+    Sigma = real(Sigma)
+    var = diag(Sigma)
+    sigma = sqrt.(var)
+    rho = Sigma ./ (sigma * sigma')
+
+    return rho, Sigma, var
+end

Computed Data/EppsCorrection/GBMwExpSaples2.jld

@@ -0,0 +1,105 @@
+## Author: Patrick Chang
+# Script file for the MM Complex Fourier Transform (Dirichlet)
+# Supporting Algorithms are at the start of the script
+#  Include:
+#           - Scale function to re-scale time to [0, 2 \pi]
+# Number of Fourier Coefficients automatically chosen so that events
+# are not aliased
+
+#---------------------------------------------------------------------------
+
+### Data Format:
+## p1;p2 = vector of stock prices - can include or not include NaNs
+## t1;t2 = vector of trading times - can include or not include NaNs
+## N = Optional input for cutoff frequency
+
+#---------------------------------------------------------------------------
+
+using ArgCheck; using LinearAlgebra
+
+#---------------------------------------------------------------------------
+### Supporting functions
+
+function scaleV2(t1, t2)
+    maxt = maximum(filter(!isnan, [t1;t2]))
+    mint = minimum(filter(!isnan, [t1;t2]))
+
+    tau1 = (2*pi) .* (t1 .- mint) ./ (maxt - mint)
+    tau2 = (2*pi) .* (t2 .- mint) ./ (maxt - mint)
+    return tau1, tau2
+end
+
+function KanataniDK(t1, t2, N)
+    t1 = filter(!isnan, t1)
+    t2 = filter(!isnan, t2)
+    L1 = length(t1)
+    L2 = length(t2)
+
+    W = zeros(L1-1, L2-1)
+
+    for i in 1:(L1-1)
+        for j in 1:(L2-1)
+            if (t1[i] == t2[j])
+                W[i,j] = 1
+            else
+                W[i,j] = (1/(2*N+1)) * sin((N + 0.5) * (t1[i] - t2[j]))/sin((t1[i] - t2[j])/2)
+            end
+        end
+    end
+    return W
+end
+
+#---------------------------------------------------------------------------
+
+# Kanatani Weight matrix implementaion using the Dirichlet Kernel   
+
+function KANcorrDK(p1, p2, t1, t2; kwargs...)
+    t1 = filter(!isnan, t1)
+    t2 = filter(!isnan, t2)
+
+    p1 = filter(!isnan, p1)
+    p2 = filter(!isnan, p2)
+
+    # Re-scale trading times
+    tau = scaleV2(t1, t2)
+    tau1 = tau[1]
+    tau2 = tau[2]
+    # Computing minimum time change
+    # minumum step size to avoid smoothing
+    dtau1 = diff(tau1)
+    dtau2 = diff(tau2)
+    taumin = minimum([dtau1; dtau2])
+    taumax = 2*pi
+    # Sampling Freq.
+    N0 = taumax/taumin
+
+    # Optional Cutoff - if not specified we use Nyquist Cutoff
+    kwargs = Dict(kwargs)
+
+    if haskey(kwargs, :N)
+        N = kwargs[:N]
+    else
+        N = trunc(Int, floor(N0/2))
+    end
+
+    #------------------------------------------------------
+
+    DiffP1 = diff(log.(p1))
+    DiffP2 = diff(log.(p2))
+
+    W1 = KanataniDK(tau1, tau1, N)
+    W2 = KanataniDK(tau2, tau2, N)
+    W12 = KanataniDK(tau1, tau2, N)
+
+    Sigma = zeros(Float64, 2, 2)
+    Sigma[1, 1] = DiffP1' * W1 * DiffP1
+    Sigma[2, 2] = DiffP2' * W2 * DiffP2
+    Sigma[1, 2] = DiffP1' * W12 * DiffP2
+    Sigma[2, 1] = Sigma[1, 2]
+
+    var = diag(Sigma)
+    sigma = sqrt.(var)
+    rho = Sigma ./ (sigma * sigma')
+
+    return rho, Sigma, var
+end