Volatility Features

fast_volatility.py

import numba
import numpy as np

@numba.njit
def sliding_window_view(a, w):
    s = a.strides[0]
    shape = a.shape[0] - w + 1, w
    return np.lib.stride_tricks.as_strided(a, shape, (s, s))

@numba.njit
# preallocate empty array and assign slice by chrisaycock
def np_shift(arr, num, fill_value=np.nan):
    result = np.empty_like(arr)
    if num > 0:
        result[:num] = fill_value
        result[num:] = arr[:-num]
    elif num < 0:
        result[num:] = fill_value
        result[:num] = arr[-num:]
    else:
        result[:] = arr
    return result

# pretty fast rolling nanstd
def nanstd(a, W):
    k = np.ones(W, dtype=int)
    m = ~np.isnan(a)
    a0 = np.where(m, a,0)

    n = np.convolve(m,k,'valid')
    c1 = np.convolve(a0, k,'valid')
    f2 = c1**2
    p2 = f2/n**2
    f1 = np.convolve((a0**2)*m,k,'valid')+n*p2
    out = np.sqrt((f1 - (2/n)*f2)/n)
    return out

# Realized Volatility: Close-Close
def realized(close: np.ndarray, window: int = 30, trading_periods: int = 365) -> np.ndarray:
    log_return = np.log(close / np_shift(close, 1))

    std = nanstd(log_return, window)
    std = np.concatenate((np.full((log_return.shape[0] - std.shape[0]), np.nan), std))
    result = std  * np.sqrt(trading_periods)
    return result

# Parkinson Volatility: High-Low Volatility
def parkinson(high: np.ndarray, low: np.ndarray, window: int = 30, trading_periods: int = 365) -> np.ndarray:
    
    rs = (1.0 / (4.0 * np.log(2.0))) * (np.log(high / low))**2.0
    
    mean = np.mean(sliding_window_view(rs, window), axis=1)
    mean = np.concatenate((np.full((high.shape[0] - mean.shape[0]), np.nan), mean))

    result = (trading_periods * mean)**0.5

    return result

# Garman-Klass Volatility: OHLC volatility
def garman_class(open: np.ndarray, high: np.ndarray, low: np.ndarray, close: np.ndarray, window: int = 30, trading_periods: int = 365) -> np.ndarray:
    
    log_hl = np.log(high / low)
    log_co = np.log(close / open)

    rs = 0.5 * log_hl**2 - (2*np.log(2)-1) * log_co**2
    
    mean = np.mean(sliding_window_view(rs, window), axis=1)
    mean = np.concatenate((np.full((open.shape[0] - mean.shape[0]), np.nan), mean))

    result = (trading_periods * mean)**0.5

    return result

# Roger-Satchell Volatility: OHLC Volatility
def rogers_satchell(open: np.ndarray, high: np.ndarray, low: np.ndarray, close: np.ndarray, window: int = 30, trading_periods: int = 365) ->  np.ndarray:
    
    log_ho = np.log(high / open)
    log_lo = np.log(low / open)
    log_co = np.log(close / open)
    
    rs = log_ho * (log_ho - log_co) + log_lo * (log_lo - log_co)
    
    mean = np.mean(sliding_window_view(rs, window), axis=1)
    mean = np.concatenate((np.full((open.shape[0] - mean.shape[0]), np.nan), mean))

    result = (trading_periods * mean)**0.5

    return result

# Garman-Klass-Yang-Zhang Volatility: OHLC Volatility
def garkla_yangzh(open : np.ndarray, high : np.ndarray, low : np.ndarray, close : np.ndarray, window: int = 30, trading_periods: int = 365) -> np.ndarray:
     
    opcl_1 = np.log(open / np_shift(close, 1))
    hilo = np.log(high / low)
    clop = np.log(close / open)

    rs = opcl_1 ** 2 + 0.5 * hilo ** 2 - (2 * np.log(2) - 1) * clop ** 2

    mean = np.mean(sliding_window_view(rs, window), axis=1)
    mean = np.concatenate((np.full((open.shape[0] - mean.shape[0]), np.nan), mean))

    result = (trading_periods * mean)**0.5

    return result

# Yang-Zhang Volatility: OHLC Volatility
def yang_zhang(open: np.ndarray, high : np.ndarray, low : np.ndarray, close : np.ndarray, window: int = 30, trading_periods: int = 365) -> np.ndarray:
     
    close_shift_1 = np_shift(close, 1)

    log_ho = np.log(high / open)
    log_lo = np.log(low / open)
    log_co = np.log(close / open)
    
    log_oc = np.log(open / close_shift_1)
    log_oc_sq = log_oc**2
    
    log_cc = np.log(close / close_shift_1)
    log_cc_sq = log_cc**2
    
    rs = log_ho * (log_ho - log_co) + log_lo * (log_lo - log_co)

    close_vol = np.sum(sliding_window_view(log_cc_sq, window), axis=1) * (1.0 / (window - 1.0))
    close_vol = np.concatenate((np.full((open.shape[0] - close_vol.shape[0]), np.nan), close_vol))
    open_vol = np.sum(sliding_window_view(log_oc_sq, window), axis=1) * (1.0 / (window - 1.0))
    open_vol = np.concatenate((np.full((open.shape[0] - open_vol.shape[0]), np.nan), open_vol))
    window_rs = np.sum(sliding_window_view(rs, window), axis=1) * (1.0 / (window - 1.0))
    window_rs = np.concatenate((np.full((open.shape[0] - window_rs.shape[0]), np.nan), window_rs))

    k = 0.34 / (1.34 + (window + 1) / (window - 1))
    result = np.sqrt(open_vol + k * close_vol + (1 - k) * window_rs) * np.sqrt(trading_periods)

    return result

### Hodges Tompkins Volatility
def hodges_timpkins(close : np.ndarray, window: int = 30, trading_periods: int = 365) -> np.ndarray:
        
    log_return = np.log(close / np_shift(close, 1))
    
    std = nanstd(log_return, window)
    std = np.concatenate((np.full((log_return.shape[0] - std.shape[0]), np.nan), std))

    vol = std  * np.sqrt(trading_periods)

    h = window
    count = log_return.size - np.count_nonzero(np.isnan(log_return))
    n = (count - h) + 1

    adj_factor = 1.0 / (1.0 - (h / n) + ((h**2 - 1) / (3 * n**2)))

    result = vol * adj_factor

    return result