In [1]:
get_ipython().ast_node_interactivity = 'all' import os import matplotlib.pyplot as plt import numpy as np np.set_printoptions(suppress=True) import matplotlib import math from PIL import Image import pandas as pd import random import tqdm from collections import defaultdict from scipy.optimize import * from concurrent.futures import ThreadPoolExecutor, ProcessPoolExecutor from queue import SimpleQueue import warnings warnings.filterwarnings('ignore')
In [2]:
targets = ( #'MATRIKS.AKSA', #'MATRIKS.AKSEN', #'MATRIKS.AKBNK', 'MATRIKS.ADEL', #'MATRIKS.CCOLA', #'MATRIKS.ALKIM', #'MATRIKS.ANHYT', #'MATRIKS.ALARK', #'MATRIKS.KAREL', #'MATRIKS.ISMEN', #'MATRIKS.ARCLK', #'MATRIKS.EREGL', #'MATRIKS.SISE', #'MATRIKS.PETKM', #'MATRIKS.SASA', #'MATRIKS.VESTL', 'MATRIKS.VESBE', #'MATRIKS.XAUTRY', #'MATRIKS.XAGTRY', #'MATRIKS.USDTRY', #'MATRIKS.EURTRY', #'MATRIKS.GBPTRY', #'MATRIKS.BIMAS', #'MATRIKS.TAVHL', 'MATRIKS.MTRKS', #'MATRIKS.MGROS', 'MATRIKS.SOKM', 'MATRIKS.TCELL', 'MATRIKS.GARAN', 'MATRIKS.ENJSA', 'MATRIKS.ASELS', 'MATRIKS.SAHOL', 'MATRIKS.TTKOM', 'MATRIKS.ENKAI', 'MATRIKS.AEFES', 'MATRIKS.GSDHO', 'MATRIKS.AKSEN', 'MATRIKS.TKFEN', 'MATRIKS.ALARK', 'MATRIKS.CIMSA', 'MATRIKS.ECILC', 'MATRIKS.SISE', 'MATRIKS.TTRAK', # test #'MATRIKS.DOHOL', #'MATRIKS.THYAO', #'MATRIKS.KCHOL', #'MATRIKS.EKGYO', # end test #'TEFAS.FYD', #'TEFAS.TCD', #'TEFAS.TCB', #'MATRIKS.SAHOL', #'MATRIKS.MAVI', #'MATRIKS.BUCIM', #'MATRIKS.GLDGR', #'MATRIKS.GWIND', #'MATRIKS.PGSUS', )
In [3]:
# Fetch data import miniclickhouse ch = miniclickhouse.connect(database="algotrade") qry = """ select * from daily_price_pivot where `date` > today() - 365 and `MATRIKS.CCOLA` is not null order by date asc """ symbol_data = defaultdict(lambda: []) for row in ch.execute(qry): for target in targets: symbol_data[target].append(row[target])
In [4]:
def start_from_zero(vals): return [x - vals[0] for x in vals]
In [5]:
def returns(values): values = list(values) prev = values[0] for val in values[1:]: yield (val - prev) / prev prev = val
In [6]:
def accumulate(values): t = 0 yield t for val in values: t += val yield t
In [7]:
_ = plt.figure(dpi=200) for target in targets: _ = plt.plot(list(accumulate(returns(symbol_data[target]))), label=target) _ = plt.legend()
Out:
In [8]:
def normalize(x): x = np.array(x) x = x + np.abs(np.min(x)) return x / np.sum(x)
In [9]:
def portfolio_to_prices(portfolio): assert len(portfolio) == len(targets) prices = [] for row in zip(*[symbol_data[target] for target in targets]): p = 0 for v, w in zip(portfolio, row): p += v * w prices.append(p) return prices
In [10]:
_ = plt.figure(dpi=200) _ = plt.plot(list(accumulate(returns(portfolio_to_prices([1 / len(targets)] * len(targets))))))
Out:
In [11]:
def random_weights(N): w = np.random.uniform(0, 1, N) return normalize(w) def random_portfolio(): return random_weights(len(targets)) _ = plt.figure(dpi=200) for _ in range(512): _ = plt.plot(list(accumulate(returns(portfolio_to_prices(random_portfolio())))), linewidth=0.2)
Out:
As you can see, there is a quite a bit of variance on the final results.
Let’s try to plot random portfolios to see if we can see the expected shape.
In [12]:
def volatility_flat(values): n = len(values) # values = np.interp(values, (np.min(values), np.max(values)), (0, 1)) ideal = np.linspace(values[0], values[-1], num=n) # return np.sum(np.abs(ideal - np.array(values))) / n return np.sum((ideal - np.array(values)) ** 2) / n
In [13]:
def hillclimb(f, start, iters, directions, lr): d = len(start) best = start for _ in range(iters): news = [np.random.normal(best, lr) for _ in range(directions)] best = min(news, key=f) return best
In [14]:
def coord_descent(fn, nparams, iters, step): x0 = [random.uniform(0, 1) for _ in range(nparams)] for iteration in range(iters): for param in range(nparams): candidates = [] for val in np.arange(0.01, 1, step): n = list(x0) n[param] = val candidates.append(np.array(n)) x0 = min(candidates, key=fn) return x0
In [15]:
def uniform_sphere(dimensions, radius): pn = np.random.normal(0, 1, dimensions) r = np.linalg.norm(pn) return radius / r * pn
In [16]:
def pattern_search(fn, nparams, ncandidates=16, decay=0.75, startradius=1): x0 = normalize([0.5 for _ in range(nparams)]) f = fn(x0) radius = startradius while radius > 0.001: candidates = [] for _ in range(ncandidates): candidate = x0 + uniform_sphere(nparams, radius) candidate = normalize(candidate) candidates.append(candidate) scores = list(map(fn, candidates)) ind = min(range(len(candidates)), key=lambda x: scores[x]) best = candidates[ind] bf = scores[ind] if bf < f: f = bf x0 = best else: radius *= decay return x0
In [17]:
target_returns = ["75.0", "85.0", "95.0", "100.0", "110.0"] target_returns = []
In [18]:
Xs = [] Ys = [] Cs = [] minvol_at_target = {} Q = SimpleQueue() def minvol_from_queue(): while True: try: rpct, std, params = Q.get(False) if "min" not in minvol_at_target: minvol_at_target["min"] = (std, params) else: mvs, mvp = minvol_at_target["min"] if std < mvs: minvol_at_target["min"] = (std, params) if rpct not in minvol_at_target: minvol_at_target[rpct] = (std, params) else: mvs, mvp = minvol_at_target[rpct] if std < mvs: minvol_at_target[rpct] = (std, params) except Exception: break def record_portfolio(params, save=True): params = normalize(params) values = portfolio_to_prices(params) rets = list(returns(values)) std = np.std(rets) # std = volatility_flat(list(accumulate(rets))) mean = np.mean(rets) pct = list(accumulate(rets))[-1] * 100 rpct = str(round(pct * 10) / 10) Q.put((rpct, std, params)) if save and pct > 0 and pct < 200 and std < 1: Xs.append(std) Ys.append(pct) Cs.append(mean / std * np.sqrt(249)) return values, std, mean, pct # Lowest volatility import cma for _ in range(4096): _ = record_portfolio(random_portfolio()) minvol_from_queue() def fitness(params): values, std, mean, pct = record_portfolio(params) return std #print("Initial pattern search") #_ = pattern_search(fitness, len(targets), 128, 0.7) print("Initial least-squares") _ = least_squares(fitness, minvol_at_target["min"][1]) print("Initial CMA") _ = cma.fmin(fitness, minvol_at_target["min"][1], 0.001, options={"verb_disp": 0}) print("Initial direct") _ = direct(fitness, [(0, 1)] * len(targets), locally_biased=False) minvol_from_queue() values, std, mean, pct = record_portfolio(minvol_at_target["min"][1]) print(f"Min volatility of initial search: {std} @ {pct:.2f}%") # Min vol at target def fitness(params): values, std, mean, pct = record_portfolio(params) pct = round(pct) return std + (pct - target) ** 2 for target in tqdm.tqdm(list(np.arange(60, 100, 1))): globals()["target"] = target # _ = minimize(fitness, [0.3] * len(targets)) _ = least_squares(fitness, [0.3] * len(targets), jac="3-point") # _ = coord_descent(fitness, len(targets), 8, 8) # _ = pattern_search(fitness, len(targets), 15, 0.8) # _ = direct(fitness, [(0, 1)] * len(targets), locally_biased=False) # _ = cma.fmin(fitness, [0.3] * len(targets), 0.001, options={"verb_disp": 0, "bounds": [0, 1]}) minvol_from_queue() for tgt in target_returns: globals()["target"] = float(tgt) print("Target:", tgt) _ = cma.fmin(fitness, [0.3] * len(targets), 0.001, options={"verb_disp": 0, "bounds": [0, 1]}) minvol_from_queue() #_ = plt.scatter(Xs, Ys, s=3, c=Cs, cmap="gray") Xs.clear() Ys.clear() Cs.clear() values, std, mean, pct = record_portfolio(minvol_at_target["min"][1]) print(f"Min volatility after search: {std} @ {pct:.2f}%")
Out:
Initial least-squares Initial CMA Initial direct Min volatility of initial search: 0.015215287173937646 @ 102.25%
Out:
100%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 40/40 [00:41<00:00, 1.04s/it]
Out:
Min volatility after search: 0.015215287173937646 @ 102.25%
Out:
In [19]:
Xs = [] Ys = [] Cs = [] for target in dict(minvol_at_target): mvs, mvp = minvol_at_target[target] values, std, mean, pct = record_portfolio(mvp, save=False) Xs.append(std) Ys.append(pct) Cs.append(mean / std * np.sqrt(249)) _ = plt.figure(dpi=200) _ = plt.scatter(Xs, Ys, s=3, c=Cs, cmap="copper") _ = plt.colorbar() values, std, mean, pct = record_portfolio(minvol_at_target["min"][1]) _ = plt.scatter(std, pct, marker="*", label=f"Min volatility ({pct:.2f}%)") #values, std, mean, pct = record_portfolio(minvol_at_target["min"][1]) #pct = int(round(pct)) #for target in target_returns: # values, std, mean, pct = record_portfolio(minvol_at_target[target][1]) # _ = plt.scatter(std, pct, marker="*", label=f"Min volatility ({target}%)") _ = plt.legend()
Out:
In [20]:
def print_weights(params): weights = [] for x, y in zip(targets, params): weights.append((x, y)) for x, y in sorted(weights, key= lambda x: x[1], reverse=True): print(f"{x} - {y * 100:.2f}") values, std, mean, pct = record_portfolio(minvol_at_target["min"][1]) pct = int(round(pct)) for target in ("min", *target_returns): print(f"{target}%") portfolio = minvol_at_target[target][1] print_weights(portfolio) print()
Out:
min% MATRIKS.VESBE - 29.51 MATRIKS.SOKM - 13.21 MATRIKS.TCELL - 9.43 MATRIKS.AKSEN - 7.33 MATRIKS.GARAN - 6.96 MATRIKS.ENKAI - 6.68 MATRIKS.ECILC - 6.54 MATRIKS.ADEL - 5.55 MATRIKS.AEFES - 4.76 MATRIKS.MTRKS - 3.73 MATRIKS.ALARK - 3.28 MATRIKS.TKFEN - 1.65 MATRIKS.CIMSA - 1.06 MATRIKS.TTRAK - 0.30 MATRIKS.SISE - 0.00 MATRIKS.ENJSA - 0.00 MATRIKS.GSDHO - 0.00 MATRIKS.SAHOL - 0.00 MATRIKS.TTKOM - 0.00 MATRIKS.ASELS - 0.00
In [21]:
_ = plt.figure(dpi=200) for target in ("min", *target_returns): values, std, mean, pct = record_portfolio(minvol_at_target[target][1]) _ = plt.plot(list(accumulate(returns(values))), label=f"{target}% return") _ = plt.legend()
Out:
In [22]:
_ = plt.figure(dpi=200) for sym, weight in zip(targets, minvol_at_target['min'][1]): _ = plt.plot([x * weight for x in accumulate(returns(symbol_data[sym]))], label=f"{sym} ({round(weight * 100)}%)") _ = plt.legend()
Out:
Allocate a given budget to the weights.
In [23]:
BUDGET = 45_000 target = 'min' weights = minvol_at_target[target][1] total = 0 for symbol, weight in zip(targets, weights): total += BUDGET * weight for symbol, weight in sorted(zip(targets, weights), key=lambda x: x[1], reverse=True): print(f"{symbol} - {weight * 100:.02f}% - {total * weight} TRY - {total * weight / symbol_data[symbol][-1]}")
Out:
MATRIKS.VESBE - 29.51% - 13281.345077031805 TRY - 1023.2161076295689 MATRIKS.SOKM - 13.21% - 5946.314963740591 TRY - 212.97689698211286 MATRIKS.TCELL - 9.43% - 4243.0226601318245 TRY - 119.45446678299055 MATRIKS.AKSEN - 7.33% - 3297.696874365242 TRY - 67.52041102303934 MATRIKS.GARAN - 6.96% - 3133.8291980823383 TRY - 104.60044052344253 MATRIKS.ENKAI - 6.68% - 3003.7622855841764 TRY - 106.14001009131366 MATRIKS.ECILC - 6.54% - 2943.619477837192 TRY - 89.20059023749066 MATRIKS.ADEL - 5.55% - 2497.9506723789423 TRY - 26.859684649235938 MATRIKS.AEFES - 4.76% - 2143.6205155356333 TRY - 31.293730153804866 MATRIKS.MTRKS - 3.73% - 1678.2534798657196 TRY - 35.5411579810614 MATRIKS.ALARK - 3.28% - 1473.765197509564 TRY - 17.79909658827976 MATRIKS.TKFEN - 1.65% - 741.2870641534146 TRY - 15.80569433162931 MATRIKS.CIMSA - 1.06% - 478.7728132076663 TRY - 4.797322777631927 MATRIKS.TTRAK - 0.30% - 136.75968483014321 TRY - 0.24539688647073968 MATRIKS.SISE - 0.00% - 1.673258129382668e-05 TRY - 4.085102854938154e-07 MATRIKS.ENJSA - 0.00% - 1.0018461888785301e-05 TRY - 3.0212490617567253e-07 MATRIKS.GSDHO - 0.00% - 8.078640634307484e-06 TRY - 1.8700557023859916e-06 MATRIKS.SAHOL - 0.00% - 4.997802344054734e-07 TRY - 1.161199429380747e-08 MATRIKS.TTKOM - 0.00% - 4.1629242424869473e-07 TRY - 1.9937376640263158e-08 MATRIKS.ASELS - 0.00% - 0.0 TRY - 0.0
In [24]:
sorted(list(minvol_at_target.keys()))
Out [24]:
['101.5', '101.9', '102.0', '102.1', '102.2', '102.3', '102.4', '102.5', '102.6', '102.7', '102.8', '102.9', '103.0', '103.1', '103.2', '103.3', '103.4', '103.5', '103.6', '103.7', '103.8', '103.9', '104.0', '104.1', '104.2', '104.3', '104.4', '104.5', '104.6', '104.7', '104.8', '104.9', '105.0', '105.1', '105.2', '105.3', '105.4', '105.5', '105.6', '105.7', '105.8', '105.9', '106.0', '106.1', '106.2', '106.3', '106.4', '106.5', '106.6', '106.7', '106.8', '106.9', '107.0', '107.1', '107.2', '107.3', '107.4', '107.5', '107.6', '107.7', '107.8', '107.9', '108.0', '108.1', '108.2', '108.3', '108.4', '108.5', '108.6', '108.7', '108.8', '108.9', '109.0', '109.1', '109.2', '109.3', '109.4', '109.5', '109.6', '109.7', '109.8', '109.9', '110.0', '110.1', '110.2', '110.3', '110.4', '110.5', '110.6', '110.7', '110.8', '110.9', '111.0', '111.1', '111.2', '111.3', '111.4', '111.5', '111.6', '111.7', '111.8', '111.9', '112.0', '112.1', '112.2', '112.3', '112.4', '112.5', '112.6', '112.7', '112.8', '112.9', '113.0', '113.1', '113.2', '113.3', '113.4', '113.5', '113.6', '113.7', '113.8', '113.9', '114.0', '114.1', '114.2', '114.3', '114.4', '114.5', '114.6', '114.7', '114.8', '114.9', '115.0', '115.1', '115.2', '115.3', '115.4', '115.5', '115.6', '115.7', '115.8', '115.9', '116.0', '116.1', '116.2', '116.3', '116.4', '116.5', '116.6', '116.7', '116.8', '116.9', '117.0', '117.1', '117.2', '117.3', '117.4', '117.5', '117.6', '117.7', '117.8', '117.9', '118.0', '118.1', '118.2', '118.3', '118.4', '118.5', '118.6', '118.7', '118.8', '118.9', '119.0', '119.1', '119.2', '119.3', '119.6', '119.7', '120.1', '120.2', '120.5', '120.6', '120.7', '120.8', '121.0', '121.1', '121.3', '121.4', '121.5', '121.6', '121.9', '123.4', '123.5', 'min']