feb2024_kelly.py

from typing import Tuple

import matplotlib
import numpy as np

matplotlib.use("TkAgg")
import pandas as pd
import random

target_annual_SR = 0.64
target_daily_SR = target_annual_SR / 16
target_annual_stdev = 0.16
target_annual_mean = target_annual_SR * target_annual_stdev

## both daily
target_stdev = target_annual_stdev / 16
target_mean = target_stdev * target_daily_SR

tail_daily_mean_range = [-0.05, 0.05]
tail_stdev = 0.03
prob_tail_range = [0.01, 0.2]

len_data_days = 2500

optimal_leverage = target_annual_mean / (target_annual_stdev)**2
## with SR 0.64 and std dev 16% the optimal leverage is 4; with these results we don't get beyond 3.6/4.4
## if you change any of the numbers above might need to experiment
plausible_leverages = np.arange(0.8*optimal_leverage,1.2*optimal_leverage, 0.01)

from syscore.interactive.progress_bar import progressBar


def results_for_one_sample():
    returns = single_bootstrap_equalised_SR()
    skew, lower_tail_ratio, upper_tail_ratio = skew_and_tail_ratio(returns)
    optimal_leverage, geometric_mean = optimal_leverage_and_geometric_mean_for_a_given_return_series(returns)
    return dict(optimal_leverage=optimal_leverage,
         max_geometric_mean=geometric_mean,
         upper_tail_ratio=upper_tail_ratio,
         lower_tail_ratio=lower_tail_ratio,
         skew=skew)

all_results = [results_for_one_sample() for __ in range(5000)]
df = pd.DataFrame(all_results)

## can now plot whatever against whatever
df.plot.scatter('skew', 'optimal_leverage')

def skew_and_tail_ratio(returns: np.array):
    as_series = pd.Series(returns)
    skew = as_series.skew()
    demeaned = as_series - as_series.mean()

    percentile_1 = demeaned.quantile(0.01)
    percentile_30 = demeaned.quantile(0.3)
    lower_ratio = percentile_1 / percentile_30
    lower_tail_ratio = lower_ratio / 4.43

    percentile_99 = demeaned.quantile(0.99)
    percentile_70 = demeaned.quantile(0.7)
    upper_ratio = percentile_99 / percentile_70
    upper_tail_ratio = upper_ratio / 4.43

    return skew, lower_tail_ratio, upper_tail_ratio

def optimal_leverage_and_geometric_mean_for_a_given_return_series(returns: np.array) -> Tuple[float, float]:
    series_of_geo_means = [annual_geometric_return_given_leverage(returns, leverage) for leverage in plausible_leverages]
    max_idx = series_of_geo_means.index(max(series_of_geo_means))

    return plausible_leverages[max_idx], series_of_geo_means[max_idx]


def single_bootstrap_equalised_SR() -> np.array:
    returns = single_bootstrap()
    return adjust_to_equalise_SR(returns)

def single_bootstrap() -> np.array:
    tail_mean_return = random.uniform(*tail_daily_mean_range)
    prob_tail = random.uniform(*prob_tail_range)

    list_of_returns = [single_day(tail_mean_return,  prob_tail) for __ in range(len_data_days)]

    return np.array(list_of_returns)

def single_day(tail_mean_return: float, prob_tail) -> float:
    outcome = random.random()
    if outcome<prob_tail:
        return random.gauss(tail_mean_return, tail_stdev)
    else:
        return random.gauss(target_mean, target_stdev)



def adjust_to_equalise_SR(returns: np.array) -> np.array:
    adj_for_stdev = returns * target_stdev / returns.std()
    adj_for_mean = adj_for_stdev + (target_mean - adj_for_stdev.mean())

    return adj_for_mean


def annual_geometric_return_given_leverage(returns: np.array, leverage: float) -> float:
    return annual_geometric_return_given_account_curve(
        pd.Series(returns*leverage)
    )

def annual_geometric_return_given_account_curve(account_curve: pd.Series) -> float:
    terminal_value = terminal_value_given_account_curve(account_curve)

    return annual_geometric_return_given_terminal_value(terminal_value, len(account_curve))

def terminal_value_given_account_curve(account_curve: pd.Series) -> float:
    one_plus = account_curve+1
    cumulative_returns = one_plus.cumprod()
    return cumulative_returns.iloc[-1]

import math

def annual_geometric_return_given_terminal_value(terminal_value:float, number_of_points: int) -> float:
    return (1+daily_geometric_return_given_terminal_value(terminal_value, number_of_points))**256-1

import numpy as np
def daily_geometric_return_given_terminal_value(terminal_value:float, number_of_points: int) -> float:
    try:
        return math.exp((1/number_of_points)*math.log(terminal_value))-1
    except:
        return np.nan