RMR generally should outperform the market as per the claims in the paper. Wondering why my implementation attached is not doing that well. Wondering if anybody would like to collaborate on improving this. In addition, how can RMR be converted into a market neutral strategy?

Clone Algorithm

4

Loading...

There was an error loading this backtest.

Backtest from
to
with
initial capital

Cumulative performance:

Algorithm
Benchmark

Custom data:

Total Returns

--

Alpha

--

Beta

--

Sharpe

--

Sortino

--

Max Drawdown

--

Benchmark Returns

--

Volatility

--

Returns | 1 Month | 3 Month | 6 Month | 12 Month |

Alpha | 1 Month | 3 Month | 6 Month | 12 Month |

Beta | 1 Month | 3 Month | 6 Month | 12 Month |

Sharpe | 1 Month | 3 Month | 6 Month | 12 Month |

Sortino | 1 Month | 3 Month | 6 Month | 12 Month |

Volatility | 1 Month | 3 Month | 6 Month | 12 Month |

Max Drawdown | 1 Month | 3 Month | 6 Month | 12 Month |

""" This is a template algorithm on Quantopian for you to adapt and fill in. """ import numpy as np import pandas as pd from scipy.stats import norm from quantopian.algorithm import attach_pipeline, pipeline_output from quantopian.pipeline import Pipeline from quantopian.pipeline.data.builtin import USEquityPricing from quantopian.pipeline.factors import AverageDollarVolume from quantopian.pipeline.filters.morningstar import Q500US from quantopian.pipeline.data.morningstar import earnings_report def initialize(context): """ Called once at the start of the algorithm. """ # Rebalance every day, n hour and m minutes after market open. schedule_function(rebalance, date_rules.every_day(), time_rules.market_open(hours = 0, minutes = 1)) # Record tracking variables at the end of each day. schedule_function(record_vars, date_rules.every_day(), time_rules.market_close()) # Create our dynamic stock selector. attach_pipeline(make_pipeline(), 'pipeline') context.output = None context.security_lists = None context.leaving_lists = None context.universe = None context.weigths = None def make_pipeline(): """ A function to create our dynamic stock selector (pipeline). Documentation on pipeline can be found here: https://www.quantopian.com/help#pipeline-title """ # Base universe set to the Q500US base_universe = Q500US() pipe = Pipeline( screen = base_universe ) return pipe def before_trading_start(context, data): """ Called every day before market open. """ context.output = pipeline_output('pipeline') # These are the securities that we are interested in trading each day. context.security_lists = context.output.index.tolist() context.leaving_lists = [stock for stock in context.portfolio.positions if stock not in context.security_lists] context.universe = context.security_lists + context.leaving_lists def rebalance(context, data): """ This rebalancing function is called according to our schedule_function settings. """ zeroOutMissingWeights = False weights = assign_weights(context, data) p = data.current(context.universe, 'price') pv = context.portfolio.portfolio_value def trade_stock(stock, tw): if not data.can_trade(stock): return oo = get_open_orders(stock) sp = p.loc[stock] if stock in p.index else 0.0 sq = 0.0 pw = 0.0 if (stock in context.portfolio.positions): sq = float(context.portfolio.positions[stock].amount) pw = sq * sp / pv ow = tw - pw oow = 0.0 for o in oo: remaining = o.amount - o.filled if np.sign(remaining) == np.sign(ow): oow += remaining * sp else: cancel_order(o) oo.remove(o) oow /= pv if np.sign(oow) != np.sign(ow) or (np.abs(oow) - np.abs(ow)) < 0.0001: for o in oo: cancel_order(o) order_target_percent(stock, tw) # For each security in our universe, order long or short positions according # to our current universe. for stock in context.security_lists: if stock in weights: trade_stock(stock, weights.loc[stock]) elif zeroOutMissingWeights: trade_stock(stock, 0.0) # Sell all previously held positions not in our new universe. for stock in context.leaving_lists: trade_stock(stock, 0.0) def record_vars(context, data): """ This function is called at the end of each day and plots certain variables. """ # Check how many long and short positions we have. longs = shorts = 0.0 for position in context.portfolio.positions.itervalues(): if position.amount > 0: longs += position.amount * position.last_sale_price if position.amount < 0: shorts += position.amount * position.last_sale_price exposure = context.portfolio.portfolio_value # Record and plot the leverage of our portfolio over time as well as the # long and short exposures. Even in minute mode, only the end-of-day # leverage is plotted. record(leverage = context.account.leverage, longs = longs / exposure, shorts = shorts / exposure) def handle_data(context, data): """ Called every minute. """ pass ####################################################################### def get_weights(context, data): pv = context.portfolio.portfolio_value w = pd.Series(data = 0.0, index = context.security_lists) holdings = context.portfolio.positions for position in holdings: w.loc[position] = holdings[position].amount * np.nan_to_num(holdings[position].last_sale_price) / pv return w def assign_weights(context, data, histLen = 30): """ Assign weights to securities that we want to order. """ p = data.history(context.universe, 'price', histLen, '1d').fillna(method = 'ffill').fillna(method = 'bfill').fillna(0.0).iloc[:-1, :] s = context.security_lists # w = get_weights(context, data) # # context.weigths = rmr_strategy(p, s, w) if context.weigths is None: context.weigths = pd.Series(1.0 / float(len(s)), index = s) context.weigths = rmr_strategy(p, s, context.weigths) return context.weigths ####################################################################### def rmr_strategy(prices, stocks, weights, eps = 1.0): """ Core of Robust Median Reviersion strategy implementation. :param prices: historical data in a form of numpy :param stocks: list of sid objects used in the algo :param weights: currently help portfolio weights :param eps: epsilon value :returns: new allocation for the portfolio securities """ numStocks = len(stocks) # update portfolio b_t = np.zeros(numStocks) x_tilde = np.zeros(numStocks) for i, stock in enumerate(stocks): b_t[i] = weights.loc[stock] if stock in weights.index else 1.0 / float(numStocks) x_tilde[i] = l1_median(prices.loc[:, stock]) / prices.ix[-1, stock] b_t_sum = np.sum(b_t) if b_t_sum != 0.0: b_t = np.divide(b_t, b_t_sum) x_bar = x_tilde.mean() # market relative deviation mark_rel_dev = x_tilde - x_bar # Calculate terms for lambda (lam) exp_return = np.dot(b_t, x_tilde) weight = eps - exp_return variability = (np.linalg.norm((mark_rel_dev))) ** 2 # test for divide-by-zero case if variability == 0.0: step_size = 0.0 else: step_size = np.fmax(0.0, weight / variability) b = b_t + step_size * mark_rel_dev b_norm = simplex_projection(b) # update portfolio new_weights = pd.Series(b_norm, index = stocks) return new_weights def simplex_projection(v, b=1): """ Projection vectors to the simplex domain. Implemented according to the paper: Efficient projections onto the l1-ball for learning in high dimensions, John Duchi, et al. ICML 2008. Implementation Time: 2011 June 17 by [email protected] AT pmail.ntu.edu.sg Optimization Problem: min_{w}\| w - v \|_{2}^{2} s.t. sum_{i=1}^{m}=z, w_{i}\geq 0 Input: A vector v \in R^{m}, and a scalar z > 0 (default=1) Output: Projection vector w :Example: >>> proj = simplex_projection([.4 ,.3, -.4, .5]) >>> print proj array([ 0.33333333, 0.23333333, 0. , 0.43333333]) >>> print proj.sum() 1.0 Original matlab implementation: John Duchi ([email protected]) Python-port: Copyright 2012 by Thomas Wiecki ([email protected]). """ v = np.asarray(v) p = len(v) # Sort v into u in descending order v = (v > 0) * v u = np.sort(v)[::-1] sv = np.cumsum(u) rho = np.where(u > (sv - b) / np.arange(1, p+1))[0][-1] theta = np.max([0, (sv[rho] - b) / (rho+1)]) w = (v - theta) w[w < 0] = 0 return w from scipy.optimize import minimize_scalar def l1_median(x): """ Computes L1 median (spatial median) using scipy.optimize.minimize_scalar :param x: a numpy 1D ndarray (vector) of values :returns: scalar estimate of L1 median of values """ x = np.asarray(x) x_min = float(np.amin(x)) x_max = float(np.amax(x)) res = minimize_scalar(dist_sum, bounds = (x_min, x_max), args = tuple(x), method='bounded') return res.x def dist_sum(m, *args): """ 1D sum of Euclidian distances :param m: scalar position :param *args: tuple of positions :returns: 1D sum of Euclidian distances """ return sum(abs(arg - m) for arg in args) #######################################################################