Back to Community
Fixed version of Ernie Chan's "Gold vs. gold-miners" stat arb

This is a fixed version of the Gold vs. gold-miners algorithm. Do not use this version, but see below for a fixed version that uses the batch_transform to compute the regression.

Clone Algorithm
705
Loading...
Backtest from to with initial capital
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
from zipline.transforms.utils import EventWindow
from zipline.transforms.stddev import MovingStandardDevWindow
from zipline.transforms.mavg import MovingAverageEventWindow
from zipline.utils.protocol_utils import ndict

import statsmodels.api as sm

class OLSWindow(EventWindow):
    def __init__(self, window_size, stock1, stock2):
        # Call the superclass constructor to set up base EventWindow
        # infrastructure. This call will create a window of window_size trading days -
        # holidays and weekends will be excluded from the tally.
        EventWindow.__init__(self, True, window_size, None)
        self.stock1 = stock1
        self.stock2 = stock2
        
        self.last_calc = None
        self.daycount = window_size
        self.refresh_period = window_size // 10
        # flag for signaling that the event window is full
        # i.e. we have window_size trading days of data
        self.full = False
        self.beta = None
        self.zscore = None
        self.spread_stddev = MovingStandardDevWindow(True, self.refresh_period, None)
        self.spread_mavg = MovingAverageEventWindow(['spread'], True, self.refresh_period, None)
        
    def handle_data(self, data):
        """
        New method to handle a data frame as sent to the algorithm's handle_data
        method.
        """
        
        price1 = data[self.stock1].price
        price2 = data[self.stock2].price
        spread = price1 - price2
        dt = data[self.stock1].datetime
        if dt < data[self.stock2].datetime:
            dt = data[self.stock2].datetime
            
        event = ndict({
                  'dt'    : data[self.stock1].datetime,
                  'price1': data[self.stock1].price,
                  'price2': data[self.stock2].price,
                  'spread': spread
                })
        
        self.update(event)
        
        if not self.beta:
            return
       
        current_date = data[self.stock1].datetime 
        self.spread = data[self.stock1].price - self.beta * data[self.stock2].price
        # at the moment, the standard deviation event window only looks at a price field,
        # while moving average will take a parameter for the field to average. So, our 
        # event fakes out the std deviation event window by repeating the spread in the price.
        spread_event = ndict({
                          'price' : self.spread,
                          'spread': self.spread,
                          'dt' : current_date
                          })
        # calculate the stddev of the spread
        self.spread_stddev.update(spread_event)
        # calculate the average of the spread
        self.spread_mavg.update(spread_event)
        
        spread_avg = self.spread_mavg.get_averages()['spread']
        spread_std = self.spread_stddev.get_stddev()
        if spread_std:
            self.zscore = (self.spread - spread_avg)/spread_std
        
    def handle_add(self, event):
        """
        This method is required by EventWindow, and is called for
        each new event added to the window. There isn't an incremental
        implementation for OLS available (clone and create one! We'll send you a t-shirt!), 
        so rather than recalculate every minute, we hold a value for refresh_period days
        and then recalculate.
        """
        if not self.last_calc:
            self.last_calc = event.dt
            return
        
        ols_age = event.dt - self.last_calc
        if ols_age.days >= self.refresh_period: 
            # EventWindow holds the current window of data in a list
            # called ticks. Here we are splitting it into two lists of
            # prices (floats) to pass to the Ordinary Least Squares (OLS)
            # implemented in statsmodels.
            
           
            p1 = [x.price1 for x in self.ticks]
            p2 = [x.price2 for x in self.ticks]
    
            model = sm.OLS(p1, p2)
            results = model.fit()
            self.beta = results.params[0]
            self.last_calc = event.dt
            
    def handle_remove(self, event):
        """
        This method is required by EventWindow, and is called whenever
        an event falls out of the active window. Because the OLS implementation
        we are using from statsmodels is not iterative, and is quite heavy computationally,
        we ignore this event and only periodically re-calcuate the OLS in handle_add (above)
        """
        # since an event is expiring, we know the window is full
        self.full = True
    

def initialize(context):
    context.gld = sid(26807)
    context.gdx = sid(32133)
    # there are 252 trading days per year on US Markets, so
    # 126 days is a 6 month window.
    # 6-month rolling ols calculation, recalculated every 30 days
    # ols is making a linear fit between price of gld and price of gdx.
    context.ols_window = OLSWindow(126, context.gld, context.gdx)
    # calculate the std dev over 
   
    
    # maximum total exposure (longs - shorts) in $US allowed
    context.max_notional = 1 * 1000 * 1000.0 #1M

def handle_data(context, data):
    context.ols_window.handle_data(data)
    if not context.ols_window.full:
        return  
    
    # calculate the current notional value of each position
    notional1 = context.portfolio.positions[context.gld].amount * data[context.gld].price
    notional2 = context.portfolio.positions[context.gdx].amount * data[context.gdx].price
    
    # if notional1 is zero, we don't have a bet on and we can buy or sell
    if notional1 == 0:
        can_buy = True
        can_sell = True
        spread_bet = 0.0
    else:
        # check that our spread bet has an absolute exposure within our max_notional limit.
        spread_bet = abs(notional1) + abs(notional2) * notional1/abs(notional1)
        can_buy = spread_bet < 0 or spread_bet < context.max_notional
        can_sell = spread_bet > 0 or spread_bet > -1 * context.max_notional
        
    
    # hit the escape hatch if we don't have enough data to do calculations.
    zscore = context.ols_window.zscore
    beta = context.ols_window.beta
    
    bet_shares = 5
    if zscore >= 2.0 and can_sell:
        # sell the spread, betting it will narrow since it is over 2 std deviations 
        # away from the average
        order(context.gld, -1 * bet_shares)
        order(context.gdx, bet_shares * beta)
            
    elif zscore <= -2.0 and can_buy:
        # buy the spread
        order(context.gld, bet_shares)
        order(context.gdx, -1 * beta * bet_shares)
        
    elif zscore <= 1.0 and zscore >= -1.0:
        reduce_position(context.portfolio, context.gld, data, bet_shares)
        reduce_position(context.portfolio, context.gdx, data, bet_shares)
            
            
def reduce_position(portfolio, stock, data, abs_quantity):
    """
    decrease exposure, regardless of position long/short.
    buy for a short position, sell for a long.
    """
    pos_amount = portfolio.positions[stock].amount
    if pos_amount > 0:
        order(stock, -1 * abs_quantity)
    elif pos_amount < 0:
        order(stock, abs_quantity)


This backtest was created using an older version of the backtester. Please re-run this backtest to see results using the latest backtester. Learn more about the recent changes.
There was a runtime error.
Disclaimer

The material on this website is provided for informational purposes only and does not constitute an offer to sell, a solicitation to buy, or a recommendation or endorsement for any security or strategy, nor does it constitute an offer to provide investment advisory services by Quantopian. In addition, the material offers no opinion with respect to the suitability of any security or specific investment. No information contained herein should be regarded as a suggestion to engage in or refrain from any investment-related course of action as none of Quantopian nor any of its affiliates is undertaking to provide investment advice, act as an adviser to any plan or entity subject to the Employee Retirement Income Security Act of 1974, as amended, individual retirement account or individual retirement annuity, or give advice in a fiduciary capacity with respect to the materials presented herein. If you are an individual retirement or other investor, contact your financial advisor or other fiduciary unrelated to Quantopian about whether any given investment idea, strategy, product or service described herein may be appropriate for your circumstances. All investments involve risk, including loss of principal. Quantopian makes no guarantees as to the accuracy or completeness of the views expressed in the website. The views are subject to change, and may have become unreliable for various reasons, including changes in market conditions or economic circumstances.

17 responses

Hello Thomas and others,

Any idea why the algorithm "stopped working" around mid-2012?

Most probably the uncertainty due to the Fiscal Cliff at the end of the year, money flowing into gold as a 'safe haven' but the basket of gold mining stocks GDX not seeing the same relationship due to market incorporating the increased capital gains tax into the stock price.

Thanks Lyle,

Here's a chart of the GLD & GDX prices:

Just eyeballing the charts, it looks like GDX just flattens out starting in early 2011, while GLD continues to rise (seems to be on the same general upward trend since early 2009). Kinda counter-intuitive...

Interesting analysis! I would be curious if a cointegration statistical test as described here would be sensitive to this and signal to e.g. stop trading this pair.

Thomas,

I looked over the paper a few days ago...not the sort of thing I'll be able to figure out without a significant time investment. Do you understand it? Do you think it could be implemented in Quantopian? I could contact one of the authors for guidance, but wouldn't even know what to ask at this point.

@Grant:

I actually meant to get back to you in the other thread but since this is more recent I'll just reply here.

The paper requires some background in Bayesian statistics and would definitely be a medium-scale project to implement in Quantopian. Luckily, someone else I sent the paper to is currently working on this. While he is using Mathematica for the first iteration it should be possible to port this to python (and by this extension Quantopian) without too much effort. So hopefully within a few weeks we'll have some good code to explore this further :).

@Grant:

I'm in the process of implementing the algorithms described in the Barber paper on Bayesian co-integration because I feel that they could be useful for various types of relative-value trading (e.g. simple pairs). This paper proposes a Bayesian approach to estimating whether or not two signals are co-integrated that avoids some of the pitfalls involved in using Dickey-Fuller methodology.

Initially I want to test the algorithm as described in the paper within Mathematica on simulated stock signals to verify both correctness of the implementation and to study how well it performs under the metrics of most interest to a trader. This investigation will likely suggest other ways to extend the algorithm after which it would make sense to translate the code into Python in order to test it against historical data on the Quantopian platform.

As @Thomas points out above, it will likely take us at least another few weeks to reach that point since there are more than a few steps involved here. You mention the possibility of contacting the authors for guidance. Perhaps they have some Matlab code that they are willing to share. If so, I wouldn't be surprised if it uses some of the routines from the toolbox that accompanies Barber's book on machine learning, which incidentally I recommend downloading from his website for additional background on this subject area.

Thanks rxs,

Here's a talk by the first author of the paper: http://techtalks.tv/talks/bayesian-conditional-cointegration/57453/. I listened to it...the approach is still kinda murky to me, but I'm getting a general sense.

Thomas & rxs,

Intuitively, one might expect an underlying causal relationship that would explain the common trending of GLD and GDX. Simplistically, I'd think that demand for ownership of physical gold (GLD) would result in interest in investing in companies that extract it from the earth (GDX). If so, might there be a detectable time lag between price changes in GLD and GDX? In other words, as a spread develops between GLD and GDX, is GLD driving GDX, or is there no detectable lag down to the minute level? Seems like a Quantopian algorithm could be concocted to have a look at the dynamics, right?

Grant

HI guys, I have played with cointegration in financial time series and it is a very interesting approach looking for mean reversion strategies...
Some comments come to my head reading this and the previous thread:

  • I read the code and I did not find any unit-root test such as ADF or PP or any stationarity test such as KPSS over the residual of the OLS regression. So, actually you do not know if these time series are actually cointegrated. I saw that you are proposing a bayesian approach to look for a cointegration relationship, but maybe a frequentist approach is easy to achieve using statsmodels (the last time I use it, they supported the adf test) because it is already implemented. Maybe, testing with adf test can be used as a filter to know when the pair is no longer cointegrated and not able to be traded.
  • Another idea would be to use a Johansen test to search for the cointegrated vectors... you have a lot of improvements over the adf test (also some problems) and a VECM that you can actually use to make forecasts (i think there is some work in statsmodels to have Johansen, but it is very preliminary).
  • There are some recent papers working in cointegration vectors as functions of time. In this way, the cointegration vector is no longer fixed, and you could catch /fit the dynamic of the stochastic trend underlying both time series, which would be very interesting from a forecasting point of view.

OK, just some ideas to discuss... I will try to get some improvements over this strategy with some additional code maybe implementing some of these ideas.
Cheers.

Damián.

Hi Damián,

Good to see you here :).

I think those are great ideas. Certainly using a prewritten test in statsmodel will be easier than the Bayesian cointegration test (not sure what the current status on that is; @rx?).
Can you post a paper on the cointegration of vectors as a function of time? Is the idea to look for cointegration between more than 2 stocks?

Thomas

Some interesting paper about time varying cointegration:

This is the first reference and is very old now:
http://www.jstor.org/discover/10.2307/3532993?uid=3737512&uid=2&uid=4&sid=21102185134967

Then, there was a "time of silence" and now you have some papers here:

There are more interesting papers, but I read them some years ago... I have to look for them again...

Damián.

Attached is a fixed version that also uses a batch_transform which makes the code a little nice.

Damian: Those are interesting references, I'll check them out as soon as I find some time.

Clone Algorithm
254
Loading...
Backtest from to with initial capital
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
import numpy as np
import statsmodels.api as sm

# there are 252 trading days per year on US Markets, so
# 126 days is a 6 month window.
# 6-month rolling ols calculation, recalculated every 30 days
# ols is making a linear fit between price of gld and price of gdx.
WINDOW_LENGTH = 126

@batch_transform(refresh_period=WINDOW_LENGTH//10, window_length=WINDOW_LENGTH)
def ols_transform(data, sid0, sid1): # receives constantly updated dataframe
    """Computes regression coefficient (slope and intercept)
    via Ordinary Least Squares between two SIDs.
    """
    p0 = data.price[sid0].values
    p1 = data.price[sid1].values
    slope = sm.OLS(p0, p1).fit().params[0]

    return slope

def initialize(context):
    context.gld = sid(26807)
    context.gdx = sid(32133)

    # maximum total exposure (longs - shorts) in $US allowed
    context.max_notional = 1 * 1000 * 1000.0 #1M

    context.spreads = []

def handle_data(context, data):
    ######################################################
    # 1. Compute regression coefficients between GLD and GX
    params = ols_transform(data, context.gld, context.gdx)
    if params is None: # Still in the warm-up period?
        return
    context.slope = params

    ######################################################
    # 2. Compute zscore of spread (remove mean and divide by std)
    zscore = compute_zscore(context, data)
    record(zscore=zscore)

    ######################################################
    # 3. Place orders
    place_orders(context, data, zscore)

def compute_zscore(context, data):
    """1. Compute the spread given slope and intercept.
       2. zscore the spread.
    """
    spread = data[context.gld].price - (context.slope * data[context.gdx].price)
    context.spreads.append(spread)
    zscore = (spread - np.mean(context.spreads[-WINDOW_LENGTH:])) / np.std(context.spreads[-WINDOW_LENGTH:])
    return zscore

def place_orders(context, data, zscore):
    """Buy spread if zscore is > 2, sell if zscore < .5.
    """
    # calculate the current notional value of each position
    notional1 = context.portfolio.positions[context.gld].amount * data[context.gld].price
    notional2 = context.portfolio.positions[context.gdx].amount * data[context.gdx].price

    # if notional1 is zero, we don't have a bet on and we can buy or sell
    if notional1 == 0:
        can_buy = True
        can_sell = True
        spread_bet = 0.0
    else:
        # check that our spread bet has an absolute exposure within our max_notional limit.
        spread_bet = abs(notional1) + abs(notional2) * notional1/abs(notional1)
        can_buy = spread_bet < 0 or spread_bet < context.max_notional
        can_sell = spread_bet > 0 or spread_bet > -1 * context.max_notional

    bet_shares = 5
    if zscore >= 2.0 and can_sell:
        # sell the spread, betting it will narrow since it is over 2 std deviations
        # away from the average
        order(context.gld, -1 * bet_shares)
        order(context.gdx, bet_shares * context.slope)

    elif zscore <= -2.0 and can_buy:
        # buy the spread
        order(context.gld, bet_shares)
        order(context.gdx, -1 * context.slope * bet_shares)

    elif zscore <= 1.0 and zscore >= -1.0:
        reduce_position(context, context.gld, data, bet_shares)
        reduce_position(context, context.gdx, data, bet_shares)


def reduce_position(context, stock, data, abs_quantity):
    """
    decrease exposure, regardless of position long/short.
    buy for a short position, sell for a long.
    """
    pos_amount = context.portfolio.positions[stock].amount
    if pos_amount > 0:
        order(stock, -1 * abs_quantity)
    elif pos_amount < 0:
        order(stock, abs_quantity)
This backtest was created using an older version of the backtester. Please re-run this backtest to see results using the latest backtester. Learn more about the recent changes.
There was a runtime error.

The code looks very clean... as soon as I come back home from SciPyCon Argentina Conference, I will play with this for sure...

Thanks for posting it!

Damián.

Ernie also mentioned half life in his book. I am trying to incorporate it as an exit strategy instead of zscore <= 1.0 and zscore >= -1.0 signal.
I added the calculation of half life to the algorithm, but everything is left intact.
The half life values look off and I'm still working on how to apply it. Perhaps someone can have a go at it.

Clone Algorithm
25
Loading...
Backtest from to with initial capital
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
import numpy as np
import statsmodels.api as sm

# there are 252 trading days per year on US Markets, so
# 126 days is a 6 month window.
# 6-month rolling ols calculation, recalculated every 30 days
# ols is making a linear fit between price of gld and price of gdx.
WINDOW_LENGTH = 126

@batch_transform(refresh_period=WINDOW_LENGTH//10, window_length=WINDOW_LENGTH)
def ols_transform(data, sid0, sid1): # receives constantly updated dataframe
    """Computes regression coefficient (slope and intercept)
    via Ordinary Least Squares between two SIDs.
    """
    p0 = data.price[sid0].values
    p1 = data.price[sid1].values
    slope = sm.OLS(p0, p1).fit().params[0]
    
    spread = p0 - (slope * p1) #Spread from the last 126 days
    
    prevz = spread[0:len(spread)-2] #Drop the latest spread
    z = spread[1:len(spread)-1] #Drop the earliest spread
    dz = z - prevz
    mu = np.mean(prevz) #Calculate the mean of spreads
    
    theta = sm.OLS(dz, prevz - mu).fit().params[0] #Apply OLS to get theta
    
    return [slope, theta]

def initialize(context):
    context.gld = sid(26807)
    context.gdx = sid(32133)

    # maximum total exposure (longs - shorts) in $US allowed
    context.max_notional = 1 * 1000 * 1000.0 #1M

    context.spreads = []
    set_slippage(slippage.VolumeShareSlippage(volume_limit=0.25,price_impact=0))
    set_commission(commission.PerShare(cost=0.05))

def handle_data(context, data):
    ######################################################
    # 1. Compute regression coefficients between GLD and GX
    params = ols_transform(data, context.gld, context.gdx)
    if params is None: # Still in the warm-up period?
        return
    context.slope = params[0]
    context.theta = params[1]
    half_life = -np.log(2)/context.theta

    ######################################################
    # 2. Compute zscore of spread (remove mean and divide by std)
    zscore = compute_zscore(context, data)
    record(zscore=zscore)

    ######################################################
    # 3. Place orders
    place_orders(context, data, zscore)

def compute_zscore(context, data):
    """1. Compute the spread given slope and intercept.
       2. zscore the spread.
    """
    spread = data[context.gld].price - (context.slope * data[context.gdx].price)
    context.spreads.append(spread)
    zscore = (spread - np.mean(context.spreads[-WINDOW_LENGTH:])) / np.std(context.spreads[-WINDOW_LENGTH:])
    return zscore

def place_orders(context, data, zscore):
    """Buy spread if zscore is > 2, sell if zscore < .5.
    """
    # calculate the current notional value of each position
    notional1 = context.portfolio.positions[context.gld].amount * data[context.gld].price
    notional2 = context.portfolio.positions[context.gdx].amount * data[context.gdx].price

    # if notional1 is zero, we don't have a bet on and we can buy or sell
    if notional1 == 0:
        can_buy = True
        can_sell = True
        spread_bet = 0.0
    else:
        # check that our spread bet has an absolute exposure within our max_notional limit.
        spread_bet = abs(notional1) + abs(notional2) * notional1/abs(notional1)
        can_buy = spread_bet < 0 or spread_bet < context.max_notional
        can_sell = spread_bet > 0 or spread_bet > -1 * context.max_notional
    
    bet_shares = 5
    if zscore >= 2.0 and can_sell:
        # sell the spread, betting it will narrow since it is over 2 std deviations
        # away from the average
        order(context.gld, -1 * bet_shares)
        order(context.gdx, bet_shares * context.slope)

    elif zscore <= -2.0 and can_buy:
        # buy the spread
        order(context.gld, bet_shares)
        order(context.gdx, -1 * context.slope * bet_shares)

    elif zscore <= 1.0 and zscore >= -1.0:
        reduce_position(context, context.gld, data, bet_shares)
        reduce_position(context, context.gdx, data, bet_shares)


def reduce_position(context, stock, data, abs_quantity):
    """
    decrease exposure, regardless of position long/short.
    buy for a short position, sell for a long.
    """
    pos_amount = context.portfolio.positions[stock].amount
    if pos_amount > 0:
        order(stock, -1 * abs_quantity)
    elif pos_amount < 0:
        order(stock, abs_quantity)
This backtest was created using an older version of the backtester. Please re-run this backtest to see results using the latest backtester. Learn more about the recent changes.
There was a runtime error.

This is a version of exactly closing your previous position. And One question, should we consider closing our position at the end of backtesting period in any case?

Clone Algorithm
19
Loading...
Backtest from to with initial capital
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
import numpy as np
import statsmodels.api as sm

# there are 252 trading days per year on US Markets, so
# 126 days is a 6 month window.
# 6-month rolling ols calculation, recalculated every 30 days
# ols is making a linear fit between price of gld and price of gdx.
WINDOW_LENGTH = 126

@batch_transform(refresh_period=WINDOW_LENGTH//10, window_length=WINDOW_LENGTH)
def ols_transform(data, sid0, sid1): # receives constantly updated dataframe
    """Computes regression coefficient (slope and intercept)
    via Ordinary Least Squares between two SIDs.
    """
    p0 = data.price[sid0].values
    p1 = data.price[sid1].values
    slope = sm.OLS(p0, p1).fit().params[0]

    return slope

def initialize(context):
    context.gld = sid(26807)
    context.gdx = sid(32133)

    # maximum total exposure (longs - shorts) in $US allowed
    context.max_notional = 1 * 1000 * 1000.0 #1M

    context.spreads = []

def handle_data(context, data):
    ######################################################
    # 1. Compute regression coefficients between GLD and GX
    params = ols_transform(data, context.gld, context.gdx)
    if params is None: # Still in the warm-up period?
        return
    context.slope = params

    ######################################################
    # 2. Compute zscore of spread (remove mean and divide by std)
    zscore = compute_zscore(context, data)
    record(zscore=zscore)

    ######################################################
    # 3. Place orders
    place_orders(context, data, zscore)

def compute_zscore(context, data):
    """1. Compute the spread given slope and intercept.
       2. zscore the spread.
    """
    spread = data[context.gld].price - (context.slope * data[context.gdx].price)
    context.spreads.append(spread)
    zscore = (spread - np.mean(context.spreads[-WINDOW_LENGTH:])) / np.std(context.spreads[-WINDOW_LENGTH:])
    return zscore

def place_orders(context, data, zscore):
    """Buy spread if zscore is > 2, sell if zscore < .5.
    """
    # calculate the current notional value of each position
    notional1 = context.portfolio.positions[context.gld].amount * data[context.gld].price
    notional2 = context.portfolio.positions[context.gdx].amount * data[context.gdx].price

    # if notional1 is zero, we don't have a bet on and we can buy or sell
    if notional1 == 0:
        can_buy = True
        can_sell = True
        spread_bet = 0.0
    else:
        # check that our spread bet has an absolute exposure within our max_notional limit.
        spread_bet = abs(notional1) + abs(notional2) * notional1/abs(notional1)
        can_buy = spread_bet < 0 or spread_bet < context.max_notional
        can_sell = spread_bet > 0 or spread_bet > -1 * context.max_notional

    bet_shares = 5
    if zscore >= 2.0 and can_sell:
        # sell the spread, betting it will narrow since it is over 2 std deviations
        # away from the average
        order(context.gld, -1 * bet_shares)
        order(context.gdx, bet_shares * context.slope)

    elif zscore <= -2.0 and can_buy:
        # buy the spread
        order(context.gld, bet_shares)
        order(context.gdx, -1 * context.slope * bet_shares)

    elif zscore <= 1.0 and zscore >= -1.0:
        order_target(context.gld,0)
        order_target(context.gdx,0)
There was a runtime error.

hacked another version with NUGT as the long ....

Clone Algorithm
73
Loading...
Backtest from to with initial capital
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
import numpy as np
import statsmodels.api as sm

# there are 252 trading days per year on US Markets, so
# 126 days is a 6 month window.
# 6-month rolling ols calculation, recalculated every 30 days
# ols is making a linear fit between price of gld and price of gdx.
WINDOW_LENGTH = 126

@batch_transform(refresh_period=WINDOW_LENGTH//10, window_length=WINDOW_LENGTH)
def ols_transform(data, sid0, sid1): # receives constantly updated dataframe
    """Computes regression coefficient (slope and intercept)
    via Ordinary Least Squares between two SIDs.
    """
    p0 = data.price[sid0].values
    p1 = data.price[sid1].values
    slope = sm.OLS(p0, p1).fit().params[0]

    return slope

def initialize(context):
    context.gld = sid(26807)
    context.gdx = sid(32133)
    context.nugt = sid(40553)
    context.dust = sid(40554)
    # maximum total exposure (longs - shorts) in $US allowed
    context.max_notional = 1 * 1000 * 1000.0 #1M

    context.spreads = []

def handle_data(context, data):
    ######################################################
    # 1. Compute regression coefficients between GLD and GX
    params = ols_transform(data, context.gld, context.gdx)
    if params is None: # Still in the warm-up period?
        return
    context.slope = params

    ######################################################
    # 2. Compute zscore of spread (remove mean and divide by std)
    zscore = compute_zscore(context, data)
    record(zscore=zscore)

    ######################################################
    # 3. Place orders
    place_orders(context, data, zscore)

def compute_zscore(context, data):
    """1. Compute the spread given slope and intercept.
       2. zscore the spread.
    """
    spread = data[context.gld].price - (context.slope * data[context.gdx].price)
    context.spreads.append(spread)
    zscore = (spread - np.mean(context.spreads[-WINDOW_LENGTH:])) / np.std(context.spreads[-WINDOW_LENGTH:])
    return zscore

def place_orders(context, data, zscore):
    """Buy spread if zscore is > 2, sell if zscore < .5.
    """
    # calculate the current notional value of each position
    notional1 = context.portfolio.positions[context.gld].amount * data[context.gld].price
    notional2 = context.portfolio.positions[context.nugt].amount * data[context.nugt].price

    # if notional1 is zero, we don't have a bet on and we can buy or sell
    if notional1 == 0:
        can_buy = True
        can_sell = True
        spread_bet = 0.0
    else:
        # check that our spread bet has an absolute exposure within our max_notional limit.
        spread_bet = abs(notional1) + abs(notional2) * notional1/abs(notional1)
        can_buy = spread_bet < 0 or spread_bet < context.max_notional
        can_sell = spread_bet > 0 or spread_bet > -1 * context.max_notional

    bet_shares = 3
    if zscore >= 2.0 and can_sell:
        # sell the spread, betting it will narrow since it is over 2 std deviations
        # away from the average
        order(context.gld, -1 * bet_shares)
        order(context.nugt, bet_shares * context.slope)

    elif zscore <= -2.0 and can_buy:
        # buy the spread
        order(context.gld, bet_shares)
        order(context.nugt, -1 * context.slope * bet_shares)

    elif zscore <= 1.0 and zscore >= -1.0:
        reduce_position(context, context.gld, data, bet_shares)
        reduce_position(context, context.nugt, data, bet_shares)


def reduce_position(context, stock, data, abs_quantity):
    """
    decrease exposure, regardless of position long/short.
    buy for a short position, sell for a long.
    """
    pos_amount = context.portfolio.positions[stock].amount
    if pos_amount > 0:
        order(stock, -1 * abs_quantity)
    elif pos_amount < 0:
        order(stock, abs_quantity)
There was a runtime error.