I was trying to something similar but my Python knowledge got on the way. Will be also interested in the algorithm when done.

Partly experimental code:

import numpy as np  
import pandas as pd  
import collections  
import copy

from sklearn import linear_model

def initialize(context):  
    set_commission(commission.PerShare(0.01))  
    set_slippage(slippage.VolumeShareSlippage(volume_limit=0.25, price_impact=0.01))  
    c = context  
    c.SPDR = {}

    c.SPDR['index'] = [symbol('SPY')]  
    c.SPDR['sectors'] = [symbol('XLY'), symbol('XLP'), symbol('XLE'), symbol('XLF'), symbol('XLV'), symbol('XLI'), symbol('XLB'), symbol('XLK'), symbol('XLU')]

    c.secs = c.SPDR['index'] + c.SPDR['sectors']  
    c.to_200 = pd.DataFrame({i : np.linspace(1, 200, 200) for i in c.secs})  
    c.to_100 = pd.DataFrame({i : np.linspace(1, 100, 100) for i in c.secs})  
    c.to_050 = pd.DataFrame({i : np.linspace(1, 50, 50) for i in c.secs})  
    c.to_030 = pd.DataFrame({i : np.linspace(1, 30, 30) for i in c.secs})  
    c.to_015 = pd.DataFrame({i : np.linspace(1, 15, 15) for i in c.secs})  
    c.to_007 = pd.DataFrame({i : np.linspace(1, 7, 7) for i in c.secs})  
    c.regr = linear_model.LinearRegression()  
    c.numsecs = float(len(c.secs))  
    c.margine = 0.30  
def handle_data(context, data):  
    c = context  
    p200 = history(200, '1d', 'price', ffill=True)  
    p100 = p200.iloc[-100:, :]  
    p050 = p200.iloc[-50:, :]  
    p030 = p200.iloc[-30:, :]  
    p015 = p200.iloc[-15:, :]  
    p007 = p200.iloc[-7:, :]  
    if len(p200.columns) != len(c.secs) or \  
    len(p100.columns) != len(c.secs) or \  
    len(p050.columns) != len(c.secs) or \  
    len(p030.columns) != len(c.secs) or \  
    len(p015.columns) != len(c.secs) or \  
    len(p007.columns) != len(c.secs) or \  
    len(data) != len(c.secs):  
        log.warn('No complete data')  
        return  
    lnp200 = np.log(p200)  
    lnp100 = lnp200.iloc[-100:, :]  
    lnp050 = lnp200.iloc[-50:, :]  
    lnp030 = lnp200.iloc[-30:, :]  
    lnp015 = lnp200.iloc[-15:, :]  
    lnp007 = lnp200.iloc[-7:, :]  

    def sgn(x):  
        return 0.0 if x is None or x == 0.0 else np.copysign(1.0, x)  
    def initiate_trade(s, w, force = False):  
        poss = float(c.portfolio.positions[s].amount)  
        oo = get_open_orders(s)  
        if force and len(oo) == 0:  
            order_target_percent(s, w)  
        elif poss == 0.0 and w != 0.0:  
            [cancel_order(i) for i in oo]  
            order_target_percent(s, w)  
        elif sgn(poss) != sgn(w):  
            [cancel_order(i) for i in oo]  
            order_target_percent(s, w)  

    c.to_200.index = lnp200.index  
    c.to_100.index = lnp100.index  
    c.to_050.index = lnp050.index  
    c.to_030.index = lnp030.index  
    c.to_015.index = lnp015.index  
    c.to_007.index = lnp007.index  
    regc200 = { i : c.regr.fit(pd.DataFrame(c.to_200[i], columns = [i]), pd.DataFrame(lnp200[i], columns = [i]), -1).coef_[-1][-1] for i in c.secs }  
    regc100 = { i : c.regr.fit(pd.DataFrame(c.to_100[i], columns = [i]), pd.DataFrame(lnp100[i], columns = [i]), -1).coef_[-1][-1] for i in c.secs }  
    regc050 = { i : c.regr.fit(pd.DataFrame(c.to_050[i], columns = [i]), pd.DataFrame(lnp050[i], columns = [i]), -1).coef_[-1][-1] for i in c.secs }  
    regc030 = { i : c.regr.fit(pd.DataFrame(c.to_030[i], columns = [i]), pd.DataFrame(lnp030[i], columns = [i]), -1).coef_[-1][-1] for i in c.secs }  
    regc015 = { i : c.regr.fit(pd.DataFrame(c.to_015[i], columns = [i]), pd.DataFrame(lnp015[i], columns = [i]), -1).coef_[-1][-1] for i in c.secs }  
    regc007 = { i : c.regr.fit(pd.DataFrame(c.to_007[i], columns = [i]), pd.DataFrame(lnp007[i], columns = [i]), -1).coef_[-1][-1] for i in c.secs }  

    dsig200 = { i : regc200[i] for i in c.secs }  
    dsig100 = { i : regc100[i] for i in c.secs }  
    dsig050 = { i : regc050[i] for i in c.secs }  
    dsig030 = { i : regc030[i] for i in c.secs }  
    dsig015 = { i : regc015[i] for i in c.secs }  
    dsig007 = { i : regc007[i] for i in c.secs }  

    sig200 = { i : regc200[i] - regc200[c.SPDR['index'][-1]] for i in c.secs }  
    sig100 = { i : regc100[i] - regc100[c.SPDR['index'][-1]] for i in c.secs }  
    sig050 = { i : regc050[i] - regc050[c.SPDR['index'][-1]] for i in c.secs }  
    sig030 = { i : regc030[i] - regc030[c.SPDR['index'][-1]] for i in c.secs }  
    sig015 = { i : regc015[i] - regc015[c.SPDR['index'][-1]] for i in c.secs }  
    sig007 = { i : regc007[i] - regc007[c.SPDR['index'][-1]] for i in c.secs }  
    sig200[c.SPDR['index'][-1]] = 0.0  
    sig100[c.SPDR['index'][-1]] = 0.0  
    sig050[c.SPDR['index'][-1]] = 0.0  
    sig030[c.SPDR['index'][-1]] = 0.0  
    sig015[c.SPDR['index'][-1]] = 0.0  
    sig007[c.SPDR['index'][-1]] = 0.0  

    dw200 = { i : sgn(dsig200[i]) / c.numsecs for i in c.secs }  
    dw100 = { i : sgn(dsig100[i]) / c.numsecs for i in c.secs }  
    dw050 = { i : sgn(dsig050[i]) / c.numsecs for i in c.secs }  
    dw030 = { i : sgn(dsig030[i]) / c.numsecs for i in c.secs }  
    dw015 = { i : sgn(dsig015[i]) / c.numsecs for i in c.secs }  
    dw007 = { i : sgn(dsig007[i]) / c.numsecs for i in c.secs }  

    w200 = { i : sgn(sig200[i]) / c.numsecs for i in c.secs }  
    w100 = { i : sgn(sig100[i]) / c.numsecs for i in c.secs }  
    w050 = { i : sgn(sig050[i]) / c.numsecs for i in c.secs }  
    w030 = { i : sgn(sig030[i]) / c.numsecs for i in c.secs }  
    w015 = { i : sgn(sig015[i]) / c.numsecs for i in c.secs }  
    w007 = { i : sgn(sig007[i]) / c.numsecs for i in c.secs }  
    w200[c.SPDR['index'][-1]] = 0.0  
    w100[c.SPDR['index'][-1]] = 0.0  
    w050[c.SPDR['index'][-1]] = 0.0  
    w030[c.SPDR['index'][-1]] = 0.0  
    w015[c.SPDR['index'][-1]] = 0.0  
    w007[c.SPDR['index'][-1]] = 0.0  
    w200[c.SPDR['index'][-1]] = -np.sum(w200.values())  
    w100[c.SPDR['index'][-1]] = -np.sum(w100.values())  
    w050[c.SPDR['index'][-1]] = -np.sum(w050.values())  
    w030[c.SPDR['index'][-1]] = -np.sum(w030.values())  
    w015[c.SPDR['index'][-1]] = -np.sum(w015.values())  
    w007[c.SPDR['index'][-1]] = -np.sum(w007.values())  
    #w200[c.SPDR['index'][-1]] += sgn(dw200[c.SPDR['index'][-1]])  
    #w100[c.SPDR['index'][-1]] += sgn(dw100[c.SPDR['index'][-1]])  
    #w050[c.SPDR['index'][-1]] += sgn(dw050[c.SPDR['index'][-1]])  
    #w030[c.SPDR['index'][-1]] += sgn(dw030[c.SPDR['index'][-1]])  
    #w015[c.SPDR['index'][-1]] += sgn(dw015[c.SPDR['index'][-1]])  
    #w007[c.SPDR['index'][-1]] += sgn(dw007[c.SPDR['index'][-1]])  
    w = { i : np.array([w007[i], w015[i] / 2., w030[i] / 4., w050[i] / 8., w100[i] / 16., w200[i] / 32.]) for i in c.secs }  
    dw = { i : np.array([dw007[i], dw015[i] / 2., dw030[i] / 4., dw050[i] / 8., dw100[i] * 16., dw200[i] / 32.]) for i in c.secs }  
    #fw = { i : dw007[i] * (1 + c.margine) if sgn(np.sum(w[i]) + np.sum(dw[i])) == dw007[i] else dw007[i] for i in c.SPDR['sectors'] }  
    #fw.update({ i : dw007[i] * (1 + c.margine) if sgn(np.sum(w[i]) + np.sum(dw[i])) == dw007[i] else dw007[i] for i in c.SPDR['index'] })  
    fw = { i : (sgn(np.sum(dw[i] + w[i])) / c.numsecs) if abs(np.sum(dw[i] + w[i])) > (1. / c.numsecs) else 0.0 for i in c.secs }  
    for i in c.secs:  
        wsum = np.sum(np.abs(fw.values()))  
        fw[i] = (fw[i] / wsum) if wsum != 0.0 else fw[i]  
    [initiate_trade(i, fw[i]) for i in c.SPDR['sectors']]  
    [initiate_trade(i, fw[i]) for i in c.SPDR['index']]

Terence,
Do you think you could get this strategy working (or close to it) with a single pair? If so, I could help you expand the model to support several pairs. You can check out this post to get an idea of how we would go about it. This kind of strategy will be much less messy to code if you take a class based approach. It will also be easier to pass different settings to different pairs and to know that everything is working properly. I know you are pretty new to programming, but Python's bread and butter is in it's object oriented features, it's never too early to introduce yourself to them.

EDIT:
At first glance I noticed two things.

1. In lines 68, 78, 82, and 100 you use 'return' to kick out of the for loop. My guess is that you just want to skip that pair and move onto the next one. You will want to change the 'return' to 'continue' if that's the case. Using return will kick you completely out of the handle_data function, if the first pair triggers that, the rest of them don't get checked. Switching to 'continue' will make the loop stop processing that pair and move onto the next one.

2. Your ols_transform is returning the slope and intercept of the OLS regression. That is probably messing up later calculations. See below for returning the slope only

# Select the slope only  
slope = sm.OLS(p0, p1).fit().params[sid2]

# return both the slope and intercept separately  
intercept, slope = sm.OLS(p0, p1).fit().params  

David

Hello David,

Thank you for your feedback! I wrote something with the class you provided. I figure my code is really messy and I'd restart it again. However, I am really stuck at line 33 and line 83. The reason why it doesn't run is because line 83 has an error called:
"Runtime exception: AttributeError: 'float' object has no attribute '_data'"

The program ran but it did not generate any trade, perhaps because there was no data. I would be grateful if you can tell me a bit on how I should resolve the problem with line 83.

Terence

Here is the code:
import math
import numpy as np
import scipy as sp
import datetime as dt
import statsmodels.tsa.stattools as ts
import statsmodels.api as sm

window_length = 30
lev = 1.2

def initialize(context):
context.std_multiplier = 1.5
set_commission(commission.PerTrade(cost=1.00))
context.pairs =[ [sid(24), sid(32146)],
[sid(46216), sid(46220)],
[sid(14517), sid(14516)],
[sid(5484), sid(6119)],
[sid(24), sid(5773)],
[sid(46170), sid(46220)],
[sid(46170), sid(46222)],
[sid(46216), sid(46222)],
[sid(46216), sid(46220)]
]

# Set the allocation per stock  
pct_per_algo = 1.0 / (2*len(context.pairs))  

# Make a separate algo for each stock.  
context.algo = [PairTradeSpread(context.pairs[t][0], context.pairs[t][1], pct_per_algo) for t in range (0, (len(context.pairs)))]

def handle_data(context, data):
for algo in context.algo:
algo.handle_data(context, data)

class PairTradeSpread(object):

# Initialize with a single stock and assign a proportion of the account.  
def __init__(self, ticker1, ticker2, allocation):  
    self.ticker1 = ticker1  
    self.ticker2 = ticker2  
    self.allocation = allocation  

@batch_transform(window_length=30, refresh_period=1)  
def avg_ratio(self, price1, price2, datapanel):  

    price1 = datapanel['price'][self.ticker1]  
    price2 = datapanel['price'][self.ticker2]  
    sid1_arr = np.array(price1)  
    sid2_arr = np.array(price2)  

    if sid1_arr is None or sid2_arr is None:  
        return None  

    return (sid1_arr/sid2_arr).mean()  

@batch_transform(window_length=30, refresh_period=1)  
def bounds(self, price1, price2, mean, context, datapanel):  

    price1 = datapanel['price'][self.ticker1]  
    price2 = datapanel['price'][self.ticker2]  
    sid1_arr = np.array(price1)  
    sid2_arr = np.array(price2)  

    if sid1_arr is None or sid2_arr is None:  
            return None  

    std_amt = (sid1_arr/sid2_arr).std()*context.std_multiplier  
    high_bound = mean + std_amt  
    low_bound = mean - std_amt  

    return (high_bound, low_bound)  

def handle_data(self, context, data):

    cash = context.portfolio.portfolio_value/(2*len(context.pairs))  
    price1 = data[self.ticker1].price  
    price2 = data[self.ticker2].price  

    current_ratio = price1/price2  
    avgratio = self.avg_ratio(price1, price2, data)  
    bounds= self.bounds(price1, price2, avgratio, context, data)  
    high_bound, low_bound = bounds  
    last_trade = None  

    if current_ratio > high_bound and last_trade is not "high":  
        order_value(self.ticker2, +lev*cash)  
        order_value(self.ticker1, -lev*cash)  
        last_trade = "high"  
    elif current_ratio < low_bound and last_trade is not "low":  
        order_value(self.ticker1, +lev*cash)  
        order_value(self.ticker2, -lev*cash)  
        last_trade = "low"  
    elif last_trade is not "mid" and last_trade is "high" and current_ratio < avgratio:  
        order_target(self.ticker1, 0)  
        order_target(self.ticker2, 0)  
        last_trade = "mid"  
    elif last_trade is not "mid" and last_trade is "low" and current_ratio > avgratio:  
        order_target(self.ticker1, 0)  
        order_target(self.ticker2, 0)  
        last_trade = "mid"  
"""  

class PairTradeZscore(object):

def __init__(self, ticker1, ticker2, allocation):  
    self.ticker1 = ticker1  
    self.ticker2 = ticker2  
    self.allocation = allocation  

def ols(self, prices, price1, price2):

    p0 = prices[price1]  
    p1 = sm.add_constant(prices[price2], prepend=True)  
    slope=sm.OLS(p0, p1).fit().params[0]  
    return slope  

def zscore(self, price1, price2, slope, context):  

    slope = self.ols(prices, price1, price2)  
    spread = price1 - (slope * price2)  
    context.spreads.append(spread)  
    zscore = (spread - np.mean(spreads[-window_length:])) / np.std(spreads[-window_length:])  
    return zscore  

def handle_data(self, context, data):  

    cash = context.portfolio.portfolio_value/(2*len(context.pairs))  
    slope = self.ols(prices, price1, price2)  
    zscore = self.zscore(price1, price2, slope, context)  

    if zscore > 2:  
        order_value(self.ticker2, +lev*cash)  
        order_value(self.ticker1, -lev*cash)  
    elif zscore < -2:  
        order_value(self.ticker1, +lev*cash)  
        order_value(self.ticker2, -lev*cash)  
    elif zscore < 0.5 and zscore > -0.5:  
        order_target(self.ticker1, 0)  
        order_target(self.ticker2, 0)  

"""  

If you format the code ```` it would be easy. Use 3 instead of 4. Ctrl + K also does the formatting.

I cleaned up your version and got it running, there is probably some conflicting orders going on because you have the same stock in several pairs. You will want to check your accounting, maybe aggregate the orders so that you can submit as few as possible.

It looks to trade every 30 min and I switched you from batch_transform back to history because batch_transform is being depreciated in favor of history. Feel free to message me if you have any questions.

David

Added cash and capital to the plots to make sure there is excess leverage. A good start I think.

Hello David and Suminda,

Thank you both for helping me out! I combined with what you both suggested and did a back-test, for a short-term period. I'm not sure whether the program will trade excessive leverage again. Please let me know how I can improve this model, I would be grateful for any feedback!

Terence

On lines 67,68, 72,73, and the other similar ones, you need to make either the numerator or denominator a float in the ratio at the end. Python 2.x uses floor division for integers, so 2/3 = 0, but 2.0/3 = 0.66666, that is going to give you odd results.

Hi David,

I edited the code and ran a lot more back-tests. I realized that the trading signal can maintain for a long period of time, making the algorithm to trade non-stop to buy over 100k, short 100k. As for line 61:

context.portfolio.cash + context.portfolio.capital_used> sum(abs(price1*position0) + abs(price2*position1) for t in range (0,)):  

I am trying to get a sum of my positions so that it does not exceed my total amount of cash and capital used, but I don't know how to define the range..

I hope you can give me some feedback on whether or not this program will still trade non-stop in the future and defining the range.

Thank you for helping me!

Sincerely,
Terence

Terence, it looks like you are getting closer. I added a few lines and record statements to help track how money is spent.

I commented out all but one pair and It looks like it only traded on that first day. I'd go through and test it with one pair while you get the PairTrade class working properly, then add the other pairs. My guess is that it goes outside of its leverage requirement and doesn't place any more trades as a result.

Hello David,

I tried again with a few more backtests. I realized that my tests would all have an error issue towards october 2006 when running the test from 2002-2014. However, there is no error issue if I run from 2006-2014. This is the backtest for one year period. It trades normally but I am sure that the leverage problem is not fixed. In this problem, I realized the leverage is only based on one pair position and its allocated cash sum. I made it so that it takes into account for the portfolio value in another version.

I wonder if you can help me on limiting gross value, based on the sum of each pair's gross position, so that the program will not trade excessively.

Thank you again for your inputs and feedback!

Sincerely,
Terence Liu

Terence,
It looks like the way the exposure and leverage are calculated take the entire portfolio into account, not just one pair. I think this might be working correctly, or really close to it, I don't see any glaring problems. It never gets more than 1.4x leveraged, which is totally reasonable, and the overall exposure to the market stays pretty low. It doesn't seem like it trades excessively either.

It looks like it's coming together well, I would say that now you should just go back and test with one pair to confirm each pair is working correctly. Then add one or two pairs at a time and make sure the transition goes as expected. If all that goes well, put any finishing touches on it and paper trade it, ideally on IB, that will give you a good idea of what to expect.

Nice work,
David

I am having a similar problem as the OP, possibly also due to thinly traded pairs. My entry/exit logic should be considerably simpler. I want to trade for reversion on any pair getting a certain distance from the prior day's closing values, closing the position if the pair gets close to reversion or close to end of day, whichever comes first.

I can successfully run the strategy with a single pair (hard-coded SIDs to context.s1 and context.s2) but whenever I run with multiple thinly traded symbols using the dict structure, I get an error that says I can't trade more than 10000000000 shares or something to that effect

Anything obvious I'm doing wrong here? What's best practice for pairs trading of multiple symbols?

Chad, it looks like you are running into thinly traded stock issues. I replaced "data[stock].close_price" with "data[stock].price" and it seems to work. When using the 'close_price' field, the prices are not forward filled, meaning that if there was no trades that minute, data[stock].close_price will be NaN. Changing to data[stock].price makes it use whatever the last trade price was.

Thanks for the fix, David. However, I tried adding some more logging and a few extra conditionals to make the exits more precise and it broke again. I did set some breakpoints and step through for a few bars and noted that the "positions" were all 0 shares, and the exception came from one of the code blocks that's supposed to flatten (order_target) to 0 shares.

In previous backtesting applications I've used, this has been the generalized structure:
1 calculate whatever signals are needed
2 Load a collection of position objects. Foreach, evaluate exit criteria and close if necessary
3 else, loop through each symbol in the universe without a position open and check entry rules

I'd love if there were a smarter way to access the collection of open positions and handle exits, generally controlling the flow smarter. Like I said, I'm new at this platform

Also, why is my return (at the point of failure) 4000%? is the market neutral position screwing up the denominator of my return calculation by making account equity equal to longPos - shortPos, and not considering the cash generated from the shorting? Small denominator would make for big RoR

Chad,
The crazy returns were due to new orders being placed when there was already open orders on the books. I made a few changes, the biggest one is that you were actually overwriting data in your main for loop. Defining "context.s1" and "context.s2" meant that the context.pClose variables were being overwritten on each pass of that for loop. Another change I made is that the algo now cancels any stale limit orders on every bar, the thought process is that it will submit a new updated order if the criteria is still met. Very cool strategy, I love that drawdown. Feel free to message me if you have further questions.

David