General Market Indicators

general_indicators.jl

# General indicators
@doc """
William’s Variable Accumulation Distribution (note may hilbert transform the price then apply this indicator)
+++ makes sense: https://www.marketvolume.com/technicalanalysis/williamsaccumulationdistribution.asp / http://www.marketinout.com/technical_analysis.php?t=Williams%27_Accumulation%7CDistribution&id=118
http://www.blastchart.com/Community/IndicatorGuide/Indicators/WilliamsAccumulationDistribution.aspx
WVAD  ( ((C-O) / (H  L)) * V )
where
C  the current period’s closing price.
O  the current period’s opening price.
H  the current period’s high price.
L  the current period’s low price.
V  the current period’s volume
https://www.incrediblecharts.com/indicators/williams_accumulation_distribution.php

https://www.incrediblecharts.com/indicators/williams_accumulate_distribute.php
- Plot 0 may be a trigger
- Sma of WVAD
- Buy or sell on zero crossings
`diffn(x::Array{Float64}; n::Int64=1)::Array{Float64}`
"""
function wvad(hlc::Matrix{Float64}; n::Int64=20)::Array{Float64}
@assert n<size(x,1) && n>0 "Argument n out of bounds."


end

@doc """
Detrended Price Oscillator
DPO = Close - Simple moving average [from (n / 2 + 1) days ago]
https://www.incrediblecharts.com/indicators/detrended_price_oscillator.php
`diffn(x::Array{Float64}; n::Int64=1)::Array{Float64}`
"""
function wvad(hlc::Matrix{Float64}; n::Int64=20)::Array{Float64}
@assert n<size(x,1) && n>0 "Argument n out of bounds."


end

@doc """
Volatility Ratio
Volatility Ratio = True Range / EMA of True Range for the past n periods
https://www.incrediblecharts.com/indicators/volatility_ratio_schwager.php
`diffn(x::Array{Float64}; n::Int64=1)::Array{Float64}`
"""
function wvad(hlc::Matrix{Float64}; n::Int64=20)::Array{Float64}
@assert n<size(x,1) && n>0 "Argument n out of bounds."


end

@doc """
Lagged differencing
`diffn(x::Array{Float64}; n::Int64=1)::Array{Float64}`
"""
function diffn(x::Array{Float64}; n::Int64=1)::Array{Float64}
    @assert n<size(x,1) && n>0 "Argument n out of bounds."
    dx = zeros(size(x))
    dx[1:n] .= NaN
    @inbounds for i=n+1:size(x,1)
        dx[i] = x[i] - x[i-n]
    end
    return dx
end

@doc """
Lagged differencing
`diffn(X::Array{Float64,2}; n::Int64=1)::Array{Float64}`
"""
function diffn(X::Array{Float64,2}; n::Int64=1)::Matrix{Float64}
    @assert n<size(X,1) && n>0 "Argument n out of bounds."
    dx = zeros(size(X))
    @inbounds for j = 1:size(X,2)
        dx[:,j] = diffn(X[:,j], n=n)
    end
    return dx
end

#=
@doc """
(Adapted from StatsBase: https://raw.githubusercontent.com/JuliaStats/StatsBase.jl/master/src/scalarstats.jl)
`Compute the mode of an arbitrary array::Array{Float64}`
"""
function mode(a::AbstractArray{Float64})
    isempty(a) && error("mode: input array cannot be empty.")
    cnts = Dict{Float64,Int64}()
    # first element
    mc = 1
    mv = a[1]
    cnts[mv] = 1
    # find the mode along with table construction
    @inbounds for i = 2 : length(a)
        x = a[i]
        if haskey(cnts, x)
            c = (cnts[x] += 1)
            if c > mc
                mc = c
                mv = x
            end
        else
            cnts[x] = 1
            # in this case: c = 1, and thus c > mc won't happen
        end
    end
    return mv
end
=#

@doc """
Compute a running or rolling arithmetic mean of an array.
`runmean(x::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}`
"""
function runmean(x::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    out[1:n-1] .= NaN
    if cumulative
        fi = 1.0:size(x,1)
        @inbounds for i = n:size(x,1)
            out[i] = sum(x[1:i])/fi[i]
        end
    else
        @inbounds for i = n:size(x,1)
            out[i] = mean(x[i-n+1:i])
        end
    end
    return out
end


@doc """
Compute a running or rolling summation of an array.
`runsum(x::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}`
"""
function runsum(x::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    if cumulative
        out = cumsum(x, dims=1)
        out[1:n-1] .= NaN
    else
        out = zeros(size(x))
        out[1:n-1] .= NaN
        @inbounds for i = n:size(x,1)
            out[i] = sum(x[i-n+1:i])
        end
    end
    return out
end

@doc """
Welles Wilder summation of an array
`wilder_sum(x::Array{Float64}; n::Int64=10)::Array{Float64}`
"""
function wilder_sum(x::Array{Float64}; n::Int64=10)::Array{Float64}
    @assert n<size(x,1) && n>0 "Argument n is out of bounds."
    nf = float(n)  # type stability -- all arithmetic done on floats
    out = zeros(size(x))
    out[1] = x[1]
    @inbounds for i = 2:size(x,1)
        out[i] = x[i] + out[i-1]*(nf-1.0)/nf
    end
    return out
end

@doc """
Compute the running or rolling mean absolute deviation of an array
`runmad(x::Array{Float64}; n::Int64=10, cumulative::Bool=true, fun::Function=median)::Array{Float64}`
"""
function runmad(x::Array{Float64}; n::Int64=10, cumulative::Bool=true, fun::Function=median)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    out[1:n-1] .= NaN
    center = 0.0
    if cumulative
        fi = collect(1.0:size(x,1))
        @inbounds for i = n:size(x,1)
            center = fun(x[1:i])
            out[i] = sum(abs.(x[1:i].-center)) / fi[i]
        end
    else
        fn = float(n)
        @inbounds for i = n:size(x,1)
            center = fun(x[i-n+1:i])
            out[i] = sum(abs.(x[i-n+1:i].-center)) / fn
        end
    end
    return out
end

@doc """
Compute the running or rolling variance of an array
`runvar(x::Array{Float64}; n::Int64=10, cumulative=true)::Array{Float64}`
"""
function runvar(x::Array{Float64}; n::Int64=10, cumulative=true)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    out[1:n-1] .= NaN
    if cumulative
        @inbounds for i = n:size(x,1)
            out[i] = var(x[1:i])
        end
    else
        @inbounds for i = n:size(x,1)
            out[i] = var(x[i-n+1:i])
        end
    end
    return out
end

@doc """
Compute the running or rolling standard deviation of an array
`runsd(x::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}`
"""
function runsd(x::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}
    return sqrt.(runvar(x, n=n, cumulative=cumulative))
end

@doc """
Compute the rolling z-score of an array
`runzscore(x::Array{Float64}; n::Int64=10)::Array{Float64}`
"""
function runzscore(x::Array{Float64,1}; n::Int64=10)#::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    m = zeros(size(x))
    sd = zeros(size(x))
    xi = zeros(size(x))
    out[1:n-1] .= NaN
        @inbounds for i = n:size(x,1)
            m[i] =  mean(x[i-n+1:i])
            sd[i] = sqrt(var(x[i-n+1:i]))
            out[i] = ((x[i] - m[i]) / sd[i])
        end
    return out
end

@doc """
Compute the running or rolling covariance of two arrays
`runcov(x::Array{Float64}, y::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}`
"""
function runcov(x::Array{Float64}, y::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}
    @assert length(x) == length(y) "Dimension mismatch: length of `x` not equal to length of `y`."
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    out[1:n-1] .= NaN
    if cumulative
        @inbounds for i = n:length(x)
            out[i] = cov(x[1:i], y[1:i])
        end
    else
        @inbounds for i = n:length(x)
            out[i] = cov(x[i-n+1:i], y[i-n+1:i])
        end
    end
    return out
end

@doc """
Compute the running or rolling correlation of two arrays
`runcor(x::Array{Float64}, y::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}`
"""
function runcor(x::Array{Float64}, y::Array{Float64}; n::Int64=10, cumulative::Bool=true)::Array{Float64}
    @assert length(x) == length(y) "Dimension mismatch: length of `x` not equal to length of `y`."
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    out[1:n-1] .= NaN
    if cumulative
        @inbounds for i = n:length(x)
            out[i] = cor(x[1:i], y[1:i])
        end
    else
        @inbounds for i = n:length(x)
            out[i] = cor(x[i-n+1:i], y[i-n+1:i])
        end
    end
    return out
end

@doc """
Compute the running or rolling maximum of an array
`runmax(x::Array{Float64}; n::Int64=10, cumulative::Bool=true, inclusive::Bool=true)::Array{Float64}`
"""
function runmax(x::Array{Float64}; n::Int64=10, cumulative::Bool=true, inclusive::Bool=true)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    if inclusive
        if cumulative
            out[n] = maximum(x[1:n])
            @inbounds for i = n+1:size(x,1)
                out[i] = max(out[i-1], x[i])
            end
        else
            @inbounds for i = n:size(x,1)
                out[i] = maximum(x[i-n+1:i])
            end
        end
        out[1:n-1] .= NaN
        return out
    else
        if cumulative
            out[n+1] = maximum(x[1:n])
            @inbounds for i = n+1:size(x,1)-1
                out[i+1] = max(out[i-1], x[i-1])
            end
        else
            @inbounds for i = n:size(x,1)-1
                out[i+1] = maximum(x[i-n+1:i])
            end
        end
        out[1:n] .= NaN
        return out
    end
end

@doc """
Compute the running or rolling minimum of an array
`runmin(x::Array{Float64}; n::Int64=10, cumulative::Bool=true, inclusive::Bool=true)::Array{Float64}`
"""
function runmin(x::Array{Float64}; n::Int64=10, cumulative::Bool=true, inclusive::Bool=true)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(size(x))
    if inclusive
        if cumulative
            out[n] = minimum(x[1:n])
            @inbounds for i = n+1:size(x,1)
                out[i] = min(out[i-1], x[i])
            end
        else
            @inbounds for i = n:size(x,1)
                out[i] = minimum(x[i-n+1:i])
            end
        end
        out[1:n-1] .= NaN
        return out
    else
        if cumulative
            out[n+1] = minimum(x[1:n])
            @inbounds for i = n+1:size(x,1)-1
                out[i+1] = min(out[i-1], x[i-1])
            end
        else
            @inbounds for i = n:size(x,1)-1
                out[i+1] = minimum(x[i-n+1:i])
            end
        end
        out[1:n] .= NaN
        return out
    end
end

@doc """
Compute the running/rolling quantile of an array
`runquantile(x::Vector{Float64}; p::Float64=0.05, n::Int=10, cumulative::Bool=true)::Vector{Float64}`
"""
function runquantile(x::Array{Float64}; p::Float64=0.05, n::Int=10, cumulative::Bool=true)::Array{Float64}
    @assert n<size(x,1) && n>1 "Argument n is out of bounds."
    out = zeros(Float64, (size(x,1), size(x,2)))
    if cumulative
        @inbounds for j in 1:size(x,2), i in 2:size(x,1)
            out[i,j] = quantile(x[1:i,j], p)
        end
        out[1,:] .= NaN
    else
        @inbounds for j in 1:size(x,2), i in n:size(x,1)
            out[i,j] = quantile(x[i-n+1:i,j], p)
        end
        out[1:n-1,:] .= NaN
    end
    return out
end


@doc """
Compute the running/rolling autocorrelation of a vector.
```
function runacf(x::Array{Float64};
                n::Int = 10,
                maxlag::Int = n-3,
                lags::AbstractArray{Int,1} = 0:maxlag,
                cumulative::Bool = true)::Matrix{Float64}
```
"""
function runacf(x::Array{Float64};
                n::Int = 10,
                maxlag::Int = n-3,
                lags::AbstractArray{Int,1} = 0:maxlag,
                cumulative::Bool = true)::Matrix{Float64}
    @assert size(x, 2) == 1 "Autocorrelation input array must be one-dimensional"
    N = size(x, 1)
    @assert n < N && n > 0
    if length(lags) == 1 && lags[1] == 0
        return ones((N, 1))
    end
    out = zeros((N, length(lags))) * NaN
    if cumulative
        @inbounds for i in n:N
            out[i,:] = acf(x[1:i], lags=lags)
        end
    else
        @inbounds for i in n:N
            out[i,:] = acf(x[i-n+1:i], lags=lags)
        end
    end
    return out
end

@doc """
Simple moving average (SMA)
`sma(x::Array{Float64}; n::Int64=10)::Array{Float64}`
"""
function sma(x::Array{Float64}; n::Int64=10)::Array{Float64}
    return runmean(x, n=n, cumulative=false)
end

@doc """
Triangular moving average (TRIMA)
`trima(x::Array{Float64}; n::Int64=10, ma::Function=sma, args...)::Array{Float64}`
"""
function trima(x::Array{Float64}; n::Int64=10, ma::Function=sma)::Array{Float64}
    return ma(ma(x, n=n), n=n)
end

@doc """
Weighted moving average (WMA)
`wma(x::Array{Float64}; n::Int64=10, wts::Array{Float64}=collect(1:n)/sum(1:n))::Array{Float64}`
"""
function wma(x::Array{Float64}; n::Int64=10, wts::Array{Float64}=collect(1:n)/sum(1:n))
    @assert n<size(x,1) && n>0 "Argument n out of bounds"
    out = fill(NaN, size(x,1))
    @inbounds for i = n:size(x,1)
        out[i] = (wts' * x[i-n+1:i])[1]
    end
    return out
end

function first_valid(x::Array{Float64})::Int
    if !isnan(x[1])
        return 1
    else
        @inbounds for i in 2:length(x)
            if !isnan(x[i])
                return i
            end
        end
    end
    return 0
end

@doc """
Exponential moving average (EMA)
`ema(x::Array{Float64}; n::Int64=10, alpha::Float64=2.0/(n+1.0), wilder::Bool=false)::Array{Float64}`
"""
function ema(x::Array{Float64}; n::Int64=10, alpha::Float64=2.0/(n+1), wilder::Bool=false)
    @assert n<size(x,1) && n>0 "Argument n out of bounds."
    if wilder
        alpha = 1.0/n
    end
    out = zeros(size(x))
    i = first_valid(x)
    out[1:n+i-2] .= NaN
    out[n+i-1] = mean(x[i:n+i-1])
    @inbounds for i = n+i:size(x,1)
        out[i] = alpha * (x[i] - out[i-1]) + out[i-1]
    end
    return out
end

@doc """
Modified moving average (MMA)
`mma(x::Array{Float64}; n::Int64=10)::Array{Float64}`
"""
function mma(x::Array{Float64}; n::Int64=10)
    return ema(x, n=n, alpha=1.0/n)
end

@doc """
Double exponential moving average (DEMA)
`dema(x::Array{Float64}; n::Int64=10, alpha=2.0/(n+1), wilder::Bool=false)::Array{Float64}`
"""
function dema(x::Array{Float64}; n::Int64=10, alpha=2.0/(n+1), wilder::Bool=false)
    return 2.0 * ema(x, n=n, alpha=alpha, wilder=wilder) -
        ema(ema(x, n=n, alpha=alpha, wilder=wilder),
            n=n, alpha=alpha, wilder=wilder)
end

@doc """
Triple exponential moving average (TEMA)
`tema(x::Array{Float64}; n::Int64=10, alpha=2.0/(n+1), wilder::Bool=false)::Array{Float64}`
"""
function tema(x::Array{Float64}; n::Int64=10, alpha=2.0/(n+1), wilder::Bool=false)
    return 3.0 * ema(x, n=n, alpha=alpha, wilder=wilder) -
        3.0 * ema(ema(x, n=n, alpha=alpha, wilder=wilder),
                  n=n, alpha=alpha, wilder=wilder) +
        ema(ema(ema(x, n=n, alpha=alpha, wilder=wilder),
                n=n, alpha=alpha, wilder=wilder),
            n=n, alpha=alpha, wilder=wilder)
end

@doc """
MESA adaptive moving average (MAMA)
`mama(x::Array{Float64}; fastlimit::Float64=0.5, slowlimit::Float64=0.05)::Matrix{Float64}`
"""
function mama(x::Array{Float64}; fastlimit::Float64=0.5, slowlimit::Float64=0.05)::Matrix{Float64}
    n = size(x,1)
    out = zeros(Float64, n, 2)
    #smooth = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    smooth_1 = 0.0
    smooth_2 = 0.0
    smooth_3 = 0.0
    smooth_4 = 0.0
    smooth_5 = 0.0
    smooth_6 = 0.0
    smooth_7 = 0.0
    #detrend = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    detrend_1 = 0.0
    detrend_2 = 0.0
    detrend_3 = 0.0
    detrend_4 = 0.0
    detrend_5 = 0.0
    detrend_6 = 0.0
    detrend_7 = 0.0
    #Q1 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    Q1_1 = 0.0
    Q1_2 = 0.0
    Q1_3 = 0.0
    Q1_4 = 0.0
    Q1_5 = 0.0
    Q1_6 = 0.0
    Q1_7 = 0.0
    #I1 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    I1_1 = 0.0
    I1_2 = 0.0
    I1_3 = 0.0
    I1_4 = 0.0
    I1_5 = 0.0
    I1_6 = 0.0
    I1_7 = 0.0
    #I2 = [0.0, 0.0]
    I2_1 = 0.0
    I2_2 = 0.0
    #Q2 = [0.0, 0.0]
    Q2_1 = 0.0
    Q2_2 = 0.0
    #Re = [0.0, 0.0]
    Re_1 = 0.0
    Re_2 = 0.0
    #Im = [0.0, 0.0]
    Im_1 = 0.0
    Im_2 = 0.0
    #per = [0.0, 0.0]
    per_1 = 0.0
    per_2 = 0.0
    #sper = [0.0, 0.0]
    sper_1 = 0.0
    sper_2 = 0.0
    #phase = [0.0, 0.0]
    phase_1 = 0.0
    phase_2 = 0.0
    jI = 0.0
    jQ = 0.0
    dphase = 0.0
    alpha = 0.0
    a = 0.0962
    b = 0.5769
    @inbounds for i = 13:n
        # Smooth and detrend price movement ====================================
        smooth_7 = (4*x[i] + 3*x[i-1] + 2*x[i-2] + x[i-3]) * 0.1
        detrend_7 = (0.0962*smooth_7+0.5769*smooth_5-0.5769*smooth_3-0.0962*smooth_1) * (0.075*per_1+0.54)
        # Compute InPhase and Quandrature components ===========================
        Q1_7 = (0.0962*detrend_7+0.5769*detrend_5-0.5769*detrend_3-0.0962*detrend_1) * (0.075*per_1+0.54)
        I1_7 = detrend_4
        # Advance phase of I1 and Q1 by 90 degrees =============================
        jQ = (0.0962*Q1_7+0.5769*Q1_5-0.5769*Q1_3-0.0962*Q1_1) * (0.075*per_1+0.54)
        jI = (0.0962*I1_7+0.5769*I1_5-0.5769*I1_3-0.0962*I1_1) * (0.075*per_1+0.54)
        # Phasor addition for 3 bar averaging ==================================
        Q2_2 = Q1_7 + jI
        I2_2 = I1_7 - jQ
        # Smooth I & Q components before applying the discriminator ============
        Q2_2 = 0.2 * Q2_2 + 0.8 * Q2_1
        I2_2 = 0.2 * I2_2 + 0.8 * I2_1
        # Homodyne discriminator ===============================================
        Re_2 = I2_2 * I2_1 + Q2_2*Q2_1
        Im_2 = I2_2 * Q2_1 - Q2_2*I2_1
        Re_2 = 0.2 * Re_2 + 0.8*Re_1
        Im_2 = 0.2 * Im_2 + 0.8*Im_1
        if (Im_2 != 0.0) & (Re_2 != 0.0)
            per_2 = 360.0/atan(Im_2/Re_2)
        end
        if per_2 > 1.5 * per_1
            per_2 = 1.5*per_1
        elseif per_2 < 0.67 * per_1
            per_2 = 0.67 * per_1
        end
        if per_2 < 6.0
            per_2 = 6.0
        elseif per_2 > 50.0
            per_2 = 50.0
        end
        per_2 = 0.2*per_2 + 0.8*per_1
        sper_2 = 0.33*per_2 + 0.67*sper_1
        if I1_7 != 0.0
            phase_2 = atan(Q1_7/I1_7)
        end
        dphase = phase_1 - phase_2
        if dphase < 1.0
            dphase = 1.0
        end
        alpha = fastlimit / dphase
        if alpha < slowlimit
            alpha = slowlimit
        end
        out[i,1] = alpha*x[i] + (1.0-alpha)*out[i-1,1]
        out[i,2] = 0.5*alpha*out[i,1] + (1.0-0.5*alpha)*out[i-1,2]
        # Reset/increment array variables
        # smooth
        smooth_1 = smooth_2
        smooth_2 = smooth_3
        smooth_3 = smooth_4
        smooth_4 = smooth_5
        smooth_5 = smooth_6
        smooth_6 = smooth_7
        # detrend
        detrend_1 = detrend_2
        detrend_2 = detrend_3
        detrend_3 = detrend_4
        detrend_4 = detrend_5
        detrend_5 = detrend_6
        detrend_6 = detrend_7
        # Q1
        Q1_1 = Q1_2
        Q1_2 = Q1_3
        Q1_3 = Q1_4
        Q1_4 = Q1_5
        Q1_5 = Q1_6
        Q1_6 = Q1_7
        # I1
        I1_1 = I1_2
        I1_2 = I1_3
        I1_3 = I1_4
        I1_4 = I1_5
        I1_5 = I1_6
        I1_6 = I1_7
        # I2
        I2_1 = I2_2
        # Q2
        Q2_1 = Q2_2
        # Re
        Re_1 = Re_2
        # Im
        Im_1 = Im_2
        # per
        per_1 = per_2
        # sper
        sper_1 = sper_2
        # phase
        phase_1 = phase_2
    end
    out[1:32,:] .= NaN
    return out
end


@doc """
Hull moving average (HMA)
`hma(x::Array{Float64}; n::Int64=20)::Array{Float64}`
"""
function hma(x::Array{Float64}; n::Int64=20)
    return wma(2 * wma(x, n=Int64(round(n/2.0))) - wma(x, n=n), n=Int64(trunc(sqrt(n))))
end

@doc """
Sine-weighted moving average
"""
function swma(x::Array{Float64}; n::Int64=10)
    @assert n<size(x,1) && n>0 "Argument n out of bounds."
    w = sin.(collect(1:n) * 180.0/6.0)  # numerator weights
    d = sum(w)  # denominator = sum(numerator weights)
    out = zeros(size(x))
    out[1:n-1] .= NaN
    @inbounds for i = n:size(x,1)
        out[i] = sum(w .* x[i-n+1:i]) / d
    end
    return out
end


@doc """
Kaufman adaptive moving average (KAMA)
"""
function kama(x::Array{Float64}; n::Int64=10, nfast::Float64=0.6667, nslow::Float64=0.0645)
    @assert n<size(x,1) && n>0 "Argument n out of bounds."
    @assert nfast>0.0 && nfast<1.0 "Argument nfast out of bounds."
    @assert nslow>0.0 && nslow<1.0 "Argument nslow out of bounds."
    dir = diffn(x, n=n)  # price direction := net change in price over past n periods
    vol = runsum(abs.(diffn(x,n=1)), n=n, cumulative=false)  # volatility/noise
    er = abs.(dir) ./ vol  # efficiency ratio
    ssc = er * (nfast-nslow) .+ nslow  # scaled smoothing constant
    sc = ssc .^ 2  # smoothing constant
    # initiliaze result variable
    out = zeros(size(x))
    i = ndims(x) > 1 ? findfirst(.!isnan.(x)).I[1] : findfirst(.!isnan.(x))
    out[1:n+i-2] .= NaN
    out[n+i-1] = mean(x[i:n+i-1])
    @inbounds for i = n+1:size(x,1)
        out[i] = out[i-1] + sc[i]*(x[i]-out[i-1])
    end
    return out
end

@doc """
Arnaud-Legoux moving average (ALMA)
`alma{Float64}(x::Array{Float64}; n::Int64=9, offset::Float64=0.85, sigma::Float64=6.0)::Array{Float64}`
"""
function alma(x::Array{Float64}; n::Int64=9, offset::Float64=0.85, sigma::Float64=6.0)
    @assert n<size(x,1) && n>0 "Argument n out of bounds."
    @assert sigma>0.0 "Argument sigma must be greater than 0."
    @assert offset>=0.0 && offset<=1 "Argument offset must be in (0,1)."
    out = zeros(size(x))
    out[1:n-1] .= NaN
    m = floor(offset*(float(n)-1.0))
    s = float(n) / sigma
    w = exp.(-(((0.0:-1.0:-float(n)+1.0) .- m) .^ 2.0) / (2.0*s*s))
    wsum = sum(w)
    if wsum != 0.0
        w = w ./ wsum
    end
    @inbounds for i = n:length(x)
        out[i] = sum(x[i-n+1] .* w)
    end
    return out
end

function lagged(x::Array{Float64}, n::Int=1)::Array{Float64}
    if n > 0
        return [fill(NaN,n); x[1:end-n]]
    elseif n < 0
        return [x[(-n+1):end]; fill(NaN,-n)]
    else
        return x
    end
end

@doc """
Zero-lag exponential moving average (ZLEMA)
`zlema(x::Array{Float64}; n::Int=10, ema_args...)::Array{Float64}`
"""
function zlema(x::Array{Float64}; n::Int=10, ema_args...)::Array{Float64}
    return ema(x+(x-lagged(x,round(Int, (n-1)/2.0))), n=n; ema_args...)
end

@doc """
Bollinger bands (moving average with standard deviation bands)
`bbands(x::Array{Float64}; n::Int64=10, sigma::Float64=2.0)::Matrix{Float64}`
*Output*
- Column 1: lower band
- Column 2: middle band
- Column 3: upper band
"""
function bbands(x::Array{Float64}; n::Int64=10, sigma::Float64=2.0, ma::Function=sma, args...)::Matrix{Float64}
    @assert n<size(x,1) && n>0 "Argument n is out of bounds."
    out = zeros(size(x,1), 3)  # cols := lower bound, ma, upper bound
    out[:,2] = ma(x, n=n, args...)
    sd = runsd(x, n=n, cumulative=false)
    out[:,1] = out[:,2] - sigma*sd
    out[:,3] = out[:,2] + sigma*sd
    return out
end

@doc """
True range
`tr(hlc::Matrix{Float64})::Array{Float64}`
"""
function tr(hlc::Matrix{Float64})::Array{Float64}
    @assert size(hlc,2) == 3 "HLC array must have 3 columns."
    n = size(hlc,1)
    out = zeros(n)
    out[1] = NaN
    @inbounds for i=2:n
        out[i] = max(hlc[i,1]-hlc[i,2], hlc[i,1]-hlc[i-1,3], hlc[i-1,3]-hlc[i,2])
    end
    return out[:,1]
end

@doc """
Average true range (uses exponential moving average)
`atr(hlc::Matrix{Float64}; n::Int64=14)::Array{Float64}`
"""
function atr(hlc::Matrix{Float64}; n::Int64=14, ma::Function=ema)::Array{Float64}
    @assert n<size(hlc,1) && n>0 "Argument n out of bounds."
    return [NaN; ma(tr(hlc)[2:end], n=n)]
end


@doc """
Keltner bands
`keltner(hlc::Matrix{Float64}; nema::Int64=20, natr::Int64=10, mult::Int64=2)::Matrix{Float64}`
*Output*
- Column 1: lower band
- Column 2: middle band
- Column 3: upper band
"""
function keltner(hlc::Array{Float64,2}; nema::Int64=20, natr::Int64=10, mult::Float64=2.0)::Matrix{Float64}
    @assert size(hlc,2) == 3 "HLC array must have 3 columns."
    out = zeros(size(hlc,1), 3)
    out[:,2] = ema(hlc[:,3], n=nema)
    out[:,1] = out[:,2] - mult*atr(hlc, n=natr)
    out[:,3] = out[:,2] + mult*atr(hlc, n=natr)
    return out
end


# Momentum-oriented technical indicator functions

@doc """
Aroon up/down/oscillator
`aroon(hl::Array{Float64,2}; n::Int64=25)::Array{Float64}`
*Output*
- Column 1: Aroon Up
- Column 2: Aroon Down
- Column 3: Aroon Oscillator
"""
function aroon(hl::Array{Float64,2}; n::Int64=25)::Matrix{Float64}
    @assert size(hl,2) == 2 "Argument `hl` must have exactly 2 columns."
    @assert n < size(hl,1) "Argument `n` must be less than the number of rows in argument `hl`."
    out = zeros(Float64, (size(hl,1),3))
    @inbounds for i in n:size(hl,1)
        out[i,1] = 100.0 * (n - findmax(hl[i-n+1:i,1])[2]) / 25.0
    end
    @inbounds for i in n:size(hl,1)
        out[i,2] = 100.0 * (n - findmin(hl[i-n+1:i,2])[2]) / 25.0
    end
    out[:,3] = out[:,1]-out[:,2]  # aroon oscillator
    out[1:n-1,:] .= NaN
    return out
end

@doc """
Donchian channel (if inclusive is set to true, will include current bar in calculations.)
`donch(hl::Array{Float64,2}; n::Int64=10, inclusive::Bool=true)::Array{Float64}`
*Output*
- Column 1: Lowest low of last `n` periods
- Column 2: Average of highest high and lowest low of last `n` periods
- Column 3: Highest high of last `n` periods
"""
function donch(hl::Array{Float64,2}; n::Int64=10, inclusive::Bool=true)::Matrix{Float64}
    @assert size(hl,2) == 2 "Argument `hl` must have exactly 2 columns."
    local lower::Array{Float64} = runmin(hl[:,2], n=n, cumulative=false, inclusive=inclusive)
    local upper::Array{Float64} = runmax(hl[:,1], n=n, cumulative=false, inclusive=inclusive)
    local middle::Array{Float64} = (lower .+ upper) ./ 2.0
    return [lower middle upper]
end

@doc """
Momentum indicator (price now vs price `n` periods back)
`momentum(x::Array{Float64}; n::Int64=1)::Array{Float64}`
"""
function momentum(x::Array{Float64}; n::Int64=1)::Array{Float64}
    @assert n>0 "Argument n must be positive."
    return diffn(x, n=n)
end

@doc """
Rate of change indicator (percent change between i'th observation and (i-n)'th observation)
`roc(x::Array{Float64}; n::Int64=1)::Array{Float64}`
"""
function roc(x::Array{Float64}; n::Int64=1)::Array{Float64}
    @assert n<size(x,1) && n>0 "Argument n out of bounds."
    out = zeros(size(x))
    @inbounds for i = n:size(x,1)
        out[i] = x[i]/x[i-n+1] - 1.0
    end
    out[1:n] .= NaN
    return out * 100.0
end

@doc """
Moving average convergence-divergence
`macd(x::Array{Float64}; nfast::Int64=12, nslow::Int64=26, nsig::Int64=9)::Array{Float64}`
*Output*
- Column 1: MACD
- Column 2: MACD Signal Line
- Column 3: MACD Histogram
"""
function macd(x::Array{Float64}; nfast::Int64=12, nslow::Int64=26, nsig::Int64=9,
              fastMA::Function=ema, slowMA::Function=ema, signalMA::Function=sma)::Matrix{Float64}
    out = zeros(Float64, (length(x),3))
    out[:,1] = fastMA(x, n=nfast) - slowMA(x, n=nslow)
    out[:,2] = signalMA(out[:,1], n=nsig)
    out[:,3] = out[:,1] - out[:,2]
    return out
end

@doc """
Relative strength index
`rsi(x::Array{Float64}; n::Int64=14, ma::Function=ema, args...)::Array{Float64}`
"""
function rsi(x::Array{Float64}; n::Int64=14, ma::Function=ema, args...)::Array{Float64}
    @assert n<size(x,1) && n>0 "Argument n is out of bounds."
    N = size(x,1)
    ups = zeros(N)
    dns = zeros(N)
    zro = 0.0
    dx = [NaN; ndims(x) > 1 ? diff(x, dims=1) : diff(x)]
    @inbounds for i=2:N
        if dx[i] > zro
            ups[i] = dx[i]
        elseif dx[i] < zro
            dns[i] = -dx[i]
        end
    end
    rs = [NaN; ma(ups[2:end], n=n; args...) ./ ma(dns[2:end], n=n; args...)]
    return 100.0 .- 100.0 ./ (1.0 .+ rs)
end

@doc """
Average directional index
`adx(hlc::Array{Float64}; n::Int64=14, wilder=true)::Array{Float64}`
*Output*
- Column 1: DI+
- Column 2: DI-
- Column 3: ADX
- Column 4: DMX
"""
# http://help.geckosoftware.com/40manual/indicators/dmi.htm
function adx(hlc::Array{Float64}; n::Int64=14, ma::Function=ema, args...)::Matrix{Float64}
    @assert n<size(hlc,1) && n>0 "Argument n is out of bounds."
    if size(hlc,2) != 3
        error("HLC array must have three columns")
    end
    N = size(hlc,1)
    updm = zeros(N)
    dndm = zeros(N)
    updm[1] = dndm[1] = NaN
    @inbounds for i = 2:N
        upmove = hlc[i,1] - hlc[i-1,1]
        dnmove = hlc[i-1,2] - hlc[i,2]
        if upmove > dnmove && upmove > 0.0
            updm[i] = upmove
        elseif dnmove > upmove && dnmove > 0.0
            dndm[i] = dnmove
        end
    end
    dip = [NaN; ma(updm[2:N], n=n; args...)] ./ atr(hlc, n=n) * 100.0 #di-
    dim = [NaN; ma(dndm[2:N], n=n; args...)] ./ atr(hlc, n=n) * 100.0 #di+
    dmx = abs.(dip-dim) ./ (dip+dim) # dmi
    adx = [fill(NaN,n); ma(dmx[n+1:N], n=n; args...)] * 100.0 # average dmi
    return [dip dim adx dmx]
end

@doc """
Parabolic stop and reverse (SAR)
`psar(hl::Array{Float64}; af_min::Float64=0.02, af_max::Float64=0.2, af_inc::Float64=af_min)::Array{Float64}`
*Arguments*
- `hl`: 2D array of high and low prices in first and second columns respectively
- `af_min`: starting/initial value for acceleration factor
- `af_max`: maximum acceleration factor (accel factor capped at this value)
- `af_inc`: increment to the acceleration factor (speed of increase in accel factor)
"""
function psar(hl::Array{Float64}; af_min::Float64=0.02, af_max::Float64=0.2, af_inc::Float64=af_min)::Array{Float64}
    @assert af_min<1.0 && af_min>0.0 "Argument af_min must be in [0,1]."
    @assert af_max<1.0 && af_max>0.0 "Argument af_max must be in [0,1]."
    @assert af_inc<1.0 && af_inc>0.0 "Argument af_inc must be in [0,1]."
    @assert size(hl,2) == 2 "Argument hl must have 2 columns."
    ls0 = 1
    ls = 0
    af0 = af_min
    af = 0.0
    ep0 = hl[1,1]
    ep = 0.0
    maxi = 0.0
    mini = 0.0
    sar = zeros(Float64,size(hl,1))
    sar[1] = hl[1,2] - std(hl[:,1]-hl[:,2])
    @inbounds for i = 2:size(hl,1)
        ls = ls0
        ep = ep0
        af = af0
        mini = min(hl[i-1,2], hl[i,2])
        maxi = max(hl[i-1,1], hl[i,1])
        # Long/short signals and local extrema
        if (ls == 1)
            ls0 = hl[i,2] > sar[i-1] ? 1 : -1
            ep0 = max(maxi, ep)
        else
            ls0 = hl[i,1] < sar[i-1] ? -1 : 1
            ep0 = min(mini, ep)
        end
        # Acceleration vector
        if ls0 == ls  # no signal change
            sar[i] = sar[i-1] + af*(ep-sar[i-1])
            af0 = (af == af_max) ? af_max : (af + af_inc)
            if ls0 == 1  # current long signal
                af0 = (ep0 > ep) ? af0 : af
                sar[i] = min(sar[i], mini)
            else  # current short signal
                af0 = (ep0 < ep) ? af0 : af
                sar[i] = max(sar[i], maxi)
            end
        else  # new signal
            af0 = af_min
            sar[i] = ep0
        end
    end
    return sar
end

@doc """
KST (Know Sure Thing) -- smoothed and summed rates of change
```
kst(x::Array{Float64};
    nroc::Array{Int64}=[10,15,20,30], navg::Array{Int64}=[10,10,10,15],
    wgts::Array{Int64}=collect(1:length(nroc)), ma::Function=sma)::Array{Float64}
```
"""
function kst(x::Array{Float64}; nroc::Array{Int64}=[10,15,20,30], navg::Array{Int64}=[10,10,10,15],
             wgts::Array{Int64}=collect(1:length(nroc)), ma::Function=sma)::Array{Float64}
    @assert length(nroc) == length(navg)
    @assert all(nroc.>0) && all(nroc.<size(x,1))
    @assert all(navg.>0) && all(navg.<size(x,1))
    N = length(x)
    k = length(nroc)
    out = zeros(size(x))
    @inbounds for j = 1:k
        out += ma(roc(x, n=nroc[j]), n=navg[j]) * wgts[j]
    end
    return out
end

@doc """
Williams %R
`wpr(hlc::Array{Float64,2}, n::Int64=14)::Array{Float64}`
"""
function wpr(hlc::Array{Float64,2}; n::Int64=14)::Array{Float64}
    hihi = runmax(hlc[:,1], n=n, cumulative=false)
    lolo = runmin(hlc[:,2], n=n, cumulative=false)
    return -100 * (hihi - hlc[:,3]) ./ (hihi - lolo)
end

@doc """
Commodity channel index
`cci(hlc::Array{Float64,2}; n::Int64=20, c::Float64=0.015, ma::Function=sma)::Array{Float64}`
"""
function cci(hlc::Array{Float64,2}; n::Int64=20, c::Float64=0.015, ma::Function=sma, args...)::Array{Float64}
    tp = (hlc[:,1] + hlc[:,2] + hlc[:,3]) ./ 3.0
    dev = runmad(tp, n=n, cumulative=false, fun=mean)
    avg = ma(tp, n=n; args...)
    return (tp - avg) ./ (c * dev)
end

@doc """
Stochastic oscillator (fast or slow)
`stoch(hlc::Array{Float64,2}; nK::Int64=14, nD::Int64=3, kind::Symbol=:fast, ma::Function=sma, args...)::Matrix{Float64}`
"""
function stoch(hlc::Array{Float64,2}; nK::Int64=14, nD::Int64=3,
               kind::Symbol=:fast, ma::Function=sma, args...)::Matrix{Float64}
    @assert kind == :fast || kind == :slow "Argument `kind` must be either :fast or :slow"
    @assert nK<size(hlc,1) && nK>0 "Argument `nK` out of bounds."
    @assert nD<size(hlc,1) && nD>0 "Argument `nD` out of bounds."
    hihi = runmax(hlc[:,1], n=nK, cumulative=false)
    lolo = runmin(hlc[:,2], n=nK, cumulative=false)
    out = zeros(Float64, (size(hlc,1),2))
    out[:,1] = (hlc[:,3]-lolo) ./ (hihi-lolo) * 100.0
    out[:,2] = ma(out[:,1], n=nD; args...)
    if kind == :slow
        out[:,1] = out[:,2]
        out[:,2] = ma(out[:,1], n=nD; args...)
    end
    return out
end

@doc """
SMI (stochastic momentum oscillator)
```
smi(hlc::Array{Float64,2}; n::Int64=13, nFast::Int64=2, nSlow::Int64=25, nSig::Int64=9,
    maFast::Function=ema, maSlow::Function=ema, maSig::Function=sma)::Matrix{Float64}
```
"""
function smi(hlc::Array{Float64,2}; n::Int64=13, nFast::Int64=2, nSlow::Int64=25, nSig::Int64=9,
             maFast::Function=ema, maSlow::Function=ema, maSig::Function=sma)::Matrix{Float64}
    hihi = runmax(hlc[:,1], n=n, cumulative=false)
    lolo = runmin(hlc[:,2], n=n, cumulative=false)
    hldif = hihi-lolo
    delta = hlc[:,3] - (hihi+lolo) / 2.0
    numer = maSlow(maFast(delta, n=nFast), n=nSlow)
    denom = maSlow(maFast(hldif, n=nFast), n=nSlow) / 2.0
    out = zeros(Float64, (size(hlc,1),2))
    out[:,1] = 100.0*(numer./denom)
    out[:,2] = maSig(out[:,1], n=nSig)
    return out
end

@doc """
Renko chart patterns
# Methods
- Traditional (Constant Box Size): `renko(x::Array{Float64}; box_size::Float64=10.0)::Array{Int}`
- ATR Dynamic Box Size: `renko(hlc::Matrix{Float64}; box_size::Float64=10.0, use_atr::Bool=false, n::Int=14)::Array{Int}`
# Output
`Array{Int}` object of size Nx1 (where N is the number rows in `x`) where each element gives the Renko bar number of the corresponding row in `x`.
"""
# Renko chart bar identification with traditional methodology (constant box size)
function renko(x::Array{Float64}; box_size::Float64=10.0)::Array{Int}
    @assert box_size != 0
    "Argument `box_size` must be nonzero."
    if box_size < 0.0
        box_size = abs(box_size)
    end
    bar_id = ones(Int, size(x,1))
    ref_pt = x[1]
    @inbounds for i in 2:size(x,1)
        if abs(x[i]-ref_pt) >= box_size
            ref_pt = x[i]
            bar_id[i:end] .+= 1
        end
    end
    return bar_id
end

# Renko chart bar identification with option to use ATR or traditional method (constant box size)
function renko(hlc::Matrix{Float64}; box_size::Float64=10.0, use_atr::Bool=false, n::Int=14)::Array{Int}
    if use_atr
        bar_id = ones(Int, size(hlc,1))
        box_sizes = atr(hlc, n=n)
        x = hlc[:,3]
        ref_pt = x[1]
        @inbounds for i in 2:size(x,1)
            if abs(x[i]-ref_pt) >= box_sizes[i]
                ref_pt = x[i]
                bar_id[i:end] .+= 1
            end
        end
        return bar_id
    else
        return renko(hlc[:,3], box_size=box_size)
    end
end


@doc """
Moving linear regression intercept (column 1) and slope (column 2)
"""
function mlr_beta(y::Array{Float64}; n::Int64=10, x::Array{Float64}=collect(1.0:n))::Matrix{Float64}
    @assert n<length(y) && n>0 "Argument n out of bounds."
    @assert size(y,2) == 1
    @assert size(x,1) == n || size(x,1) == size(y,1)
    const_x = size(x,1) == n
    out = zeros(Float64, (length(y),2))
    out[1:n-1,:] .= NaN
    xbar = mean(x)
    ybar = runmean(y, n=n, cumulative=false)
    @inbounds for i = n:length(y)
        yi = y[i-n+1:i]
        xi = const_x ? x : x[i-n+1:i]
        out[i,2] = cov(xi,yi) / var(xi)
        out[i,1] = ybar[i] - out[i,2]*xbar
    end
    return out
end

@doc """
Moving linear regression slope
"""
function mlr_slope(y::Array{Float64}; n::Int64=10, x::Array{Float64}=collect(1.0:n))::Array{Float64}
    @assert n<length(y) && n>0 "Argument n out of bounds."
    @assert size(y,2) == 1
    @assert size(x,1) == n || size(x,1) == size(y,1)
    const_x = size(x,1) == n
    out = zeros(size(y))
    out[1:n-1] .= NaN
    @inbounds for i = n:length(y)
        yi = y[i-n+1:i]
        xi = const_x ? x : x[i-n+1:i]
        out[i] = cov(xi,yi) / var(xi)
    end
    return out
end

@doc """
Moving linear regression y-intercept
"""
function mlr_intercept(y::Array{Float64}; n::Int64=10, x::Array{Float64}=collect(1.0:n))::Array{Float64}
    @assert n<length(y) && n>0 "Argument n out of bounds."
    @assert size(y,2) == 1
    @assert size(x,1) == n || size(x,1) == size(y,1)
    const_x = size(x,1) == n
    out = zeros(size(y))
    out[1:n-1] .= NaN
    xbar = mean(x)
    ybar = runmean(y, n=n, cumulative=false)
    @inbounds for i = n:length(y)
        yi = y[i-n+1:i]
        xi = const_x ? x : x[i-n+1:i]
        out[i] = ybar[i] - xbar*(cov(xi,yi)/var(xi))
    end
    return out
end

@doc """
Moving linear regression predictions
"""
function mlr(y::Array{Float64}; n::Int64=10)::Array{Float64}
    b = mlr_beta(y, n=n)
    return b[:,1] + b[:,2]*float(n)
end

@doc """
Moving linear regression standard errors
"""
function mlr_se(y::Array{Float64}; n::Int64=10)::Array{Float64}
    yhat = mlr(y, n=n)
    r = zeros(Float64, n)
    out = zeros(size(y))
    out[1:n-1] .= NaN
    nf = float(n)
    @inbounds for i = n:length(y)
        r = y[i-n+1:i] .- yhat[i]
        out[i] = sqrt(sum(r.^2)/nf)
    end
    return out
end

@doc """
Moving linear regression upper bound
"""
function mlr_ub(y::Array{Float64}; n::Int64=10, se::Float64=2.0)::Array{Float64}
    return y + se*mlr_se(y, n=n)
end

@doc """
Moving linear regression lower bound
"""
function mlr_lb(y::Array{Float64}; n::Int64=10, se::Float64=2.0)::Array{Float64}
    return y - se*mlr_se(y, n=n)
end

@doc """
Moving linear regression bands
*Output:*
Column 1: Lower bound
Column 2: Regression estimate
Column 3: Upper bound
"""
function mlr_bands(y::Array{Float64}; n::Int64=10, se::Float64=2.0)::Matrix{Float64}
    out = zeros(Float64, (length(y),3))
    out[1:n-1,:] .= NaN
    out[:,2] = mlr(y, n=n)
    out[:,1] = mlr_lb(y, n=n, se=se)
    out[:,3] = mlr_ub(y, n=n, se=se)
    return out
end

@doc """
Moving linear regression R-squared or adjusted R-squared
"""
function mlr_rsq(y::Array{Float64}; n::Int64=10, adjusted::Bool=false)::Array{Float64}
    yhat = mlr(y, n=n)
    rsq = runcor(y, yhat, n=n, cumulative=false) .^ 2
    if adjusted
        return rsq .- (1.0.-rsq)*(1.0/(float(n).-2.0))
    else
        return rsq
    end
end



# Miscellaneous utilities

@doc """
Find where `x` crosses over `y` (returns boolean vector where crossover occurs)
"""
function crossover(x::Array{Float64}, y::Array{Float64})
    @assert size(x,1) == size(y,1)
    out = falses(size(x))
    @inbounds for i in 2:size(x,1)
        out[i] = ((x[i] > y[i]) && (x[i-1] < y[i-1]))
    end
    return out
end

@doc """
Find where `x` crosses under `y` (returns boolean vector where crossunder occurs)
"""
function crossunder(x::Array{Float64}, y::Array{Float64})
    @assert size(x,1) == size(y,1)
    out = falses(size(x))
    @inbounds for i in 2:size(x,1)
        out[i] = ((x[i] < y[i]) && (x[i-1] > y[i-1]))
    end
    return out
end