Stop Loss

Today I want to share and present an example of the flexible Stop Loss functionality that I added to the Systematic Investor Toolbox.

Let’s examine a simple Moving Average Crossover strategy:

• Buy is triggered once fast moving average crosses above the slow moving average
• Sell is triggered once fast moving average crosses below the slow moving average

Here is an example using 20 and 50 day moving averages on SPY.

###############################################################################
# Load Systematic Investor Toolbox (SIT)
# https://systematicinvestor.wordpress.com/systematic-investor-toolbox/
###############################################################################
setInternet2(TRUE)
con = gzcon(url('http://www.systematicportfolio.com/sit.gz', 'rb'))
    source(con)
close(con)
 
    #*****************************************************************
    # Load historical data
    #******************************************************************  
    load.packages('quantmod')  

    tickers = spl('SPY')	
	
    data <- new.env()
        getSymbols(tickers, src = 'yahoo', from = '1970-01-01', env = data, auto.assign = T)
        for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)		
    bt.prep(data, align='keep.all', dates='1999::')

    #*****************************************************************
    # Code Strategies
    #****************************************************************** 
    prices = data$prices   
	
    models = list()
	
    #*****************************************************************
    # Code Strategies
    #****************************************************************** 
    data$weight[] = NA
        data$weight[] = 1
    models$buy.hold = bt.run.share(data, clean.signal=T)

    #*****************************************************************
    # Code Strategies : MA Cross Over
    #****************************************************************** 
    sma.fast = SMA(prices, 20)
    sma.slow = SMA(prices, 50)

    buy.signal = iif(cross.up(sma.fast, sma.slow), 1, NA)
	
    data$weight[] = NA
        data$weight[] = iif(cross.up(sma.fast, sma.slow), 1, iif(cross.dn(sma.fast, sma.slow), 0, NA))
    models$ma.cross = bt.run.share(data, clean.signal=T, trade.summary = TRUE)

Next, I created the custom.stop.fn() function in bt.stop.r at github to handle user defined stop functions. Let’s have a look at 3 different flavors of stops below:

    #*****************************************************************
    # User Defined Stops
    #****************************************************************** 
    # fixed stop: exit trade once price falls below % from entry price
    fixed.stop <- function(weight, price, tstart, tend, pstop) {
        index = tstart : tend
        if(weight > 0)
            price[ index ] < (1 - pstop) * price[ tstart ]
        else
            price[ index ] > (1 + pstop) * price[ tstart ]
    }
	
    # trailing stop: exit trade once price falls below % from max price since start of trade
    trailing.stop <- function(weight, price, tstart, tend, pstop) {
        index = tstart : tend
        if(weight > 0) 
            price[ index ] < (1 - pstop) * cummax(price[ index ])
        else
            price[ index ] > (1 + pstop) * cummin(price[ index ])
    }	
	
    # trailing stop: exit trade once price either
    # - falls below % from max price since start of trade OR
    # - rises above % from entry price
    trailing.stop.profit.target <- function(weight, price, tstart, tend, pstop, pprofit) {
        index = tstart : tend
        if(weight > 0) {
            temp = price[ index ] < (1 - pstop) * cummax(price[ index ])
			
            # profit target
            temp = temp | price[ index ] > (1 + pprofit) * price[ tstart ]
        } else {
            temp = price[ index ] > (1 + pstop) * cummin(price[ index ])
			
            # profit target
            temp = temp | price[ index ] < (1 - pprofit) * price[ tstart ]		
        }
        return( temp )	
    }

Each user defined stop function have 4 required inputs:

• weight – signal or weight indicating position
• price – price for asset
• tstart – index(time) when trade was open
• tend – index(time) when trade is closed

Additionally, you can provide any other required information.

For example, the fixed.stop function above, will exit the trade once price falls below given % (pstop) from entry price or trade is closed. The trailing.stop.profit.target function will exit the trade once either:

• trailing stop is reached. i.e. price falls below given % (pstop) from max price since start of trade OR
• profit target is reached. i.e. price rises above given % (pprofit) from entry price

Following is an example of integrating these stop losses into the back-tests. Please note that I only use the Buy rule from the original strategy to open the trades and to close the trades I use the stop loss logic.

    #*****************************************************************
    # Exit using fixed stop
    #****************************************************************** 
    data$weight[] = NA
        data$weight[] = custom.stop.fn(coredata(buy.signal), coredata(prices), fixed.stop, 
		pstop = 1/100)
    models$ma.cross.fixed.stop = bt.run.share(data, clean.signal=T, trade.summary = TRUE)
		
    #*****************************************************************
    # Exit using trailing stop
    #****************************************************************** 
    data$weight[] = NA
    	data$weight[] = custom.stop.fn(coredata(buy.signal), coredata(prices), trailing.stop, 
		pstop = 1/100)
    models$ma.cross.trailing.stop = bt.run.share(data, clean.signal=T, trade.summary = TRUE)
	
    #*****************************************************************
    # Exit using trailing stop or profit target
    #****************************************************************** 		
    data$weight[] = NA
    	data$weight[] = custom.stop.fn(coredata(buy.signal), coredata(prices), trailing.stop.profit.target, 
		pstop = 1/100, pprofit = 1.5/100)
    models$ma.cross.trailing.stop.profit.target = bt.run.share(data, clean.signal=T, trade.summary = TRUE)
	
    #*****************************************************************
    # Create Report
    #****************************************************************** 	
    strategy.performance.snapshoot(models, T)

Nothing exciting to report. The fixed stop loss is not triggered often and resembles the Buy and Hold. The trailing stops reduce the time strategy is in the market and also reduce returns. The trailing stop with profit target is doing a bit better.

We can see functionality of the stops more clearly in the picture below:

    #*****************************************************************
    # Create Plot
    #****************************************************************** 	
    dates = '2010::2010'
    # add moving averages to the strategy plot
    extra.plot.fn <- function() {
	plota.lines(sma.fast, col='red')
	plota.lines(sma.slow, col='blue')
    }
	
    layout(1:4)
    bt.stop.strategy.plot(data, models$ma.cross, dates = dates, layout=T, main = 'MA Cross', extra.plot.fn = extra.plot.fn, plotX = F)
    bt.stop.strategy.plot(data, models$ma.cross.fixed.stop, dates = dates, layout=T, main = 'Fixed Stop', plotX = F)
    bt.stop.strategy.plot(data, models$ma.cross.trailing.stop, dates = dates, layout=T, main = 'Trailing Stop', plotX = F)
    bt.stop.strategy.plot(data, models$ma.cross.trailing.stop.profit.target, dates = dates, layout=T, main = 'Trailing Stop and Profit Target')	

The shaded green area indicates the times when strategy is invested. In the top chart, the strategy is long once fast moving average (red line) is above the slow moving average (blue line).

The next chart (Fixed Stop) is fully invested because market has taken off since we entered the trade and stop, which is % relative to the trade entry price, is never triggered.

The last 2 charts show the trailing stop logic. Strategy is less often invested once we enforce the aggressive profit target.

Have fun with creating your own customized stop loss functions and let me know if you run into problems.

To view the complete source code for this example, please have a look at the bt.stop.ma.cross.test() function in bt.stop.test.r at github.

Advertisement

Stochastic Oscillator

I came across the link to the John Ehlers paper: Predictive Indicators for Effective Trading Strategies, while reading the Dekalog Blog. John Ehlers offers a different way to smooth prices and incorporate the new filter into the oscillator construction. Fortunately, the EasyLanguage code was also provided and i was able to translate it into R.

Following are some results from the paper and test of John’s stochastic oscillator.

############################################################################### 
# Super Smoother filter 2013 John F. Ehlers
# http://www.stockspotter.com/files/PredictiveIndicators.pdf
############################################################################### 
super.smoother.filter <- function(x) {
    a1 = exp( -1.414 * pi / 10.0 )
    b1 = 2.0 * a1 * cos( (1.414*180.0/10.0) * pi / 180.0 )
    c2 = b1
    c3 = -a1 * a1
    c1 = 1.0 - c2 - c3

    x = c1 * (x + mlag(x)) / 2
        x[1] = x[2]
 
    out = x * NA
        out[] = filter(x, c(c2, c3), method='recursive', init=c(0,0))
    out
}

# Roofing filter 2013 John F. Ehlers
roofing.filter <- function(x) {
    # Highpass filter cyclic components whose periods are shorter than 48 bars
    alpha1 = (cos((0.707*360 / 48) * pi / 180.0 ) + sin((0.707*360 / 48) * pi / 180.0 ) - 1) / cos((0.707*360 / 48) * pi / 180.0 )
    
    x = (1 - alpha1 / 2)*(1 - alpha1 / 2)*( x - 2*mlag(x) + mlag(x,2))
        x[1] = x[2] = x[3]
    
    HP = x * NA
    HP[] = filter(x, c(2*(1 - alpha1), - (1 - alpha1)*(1 - alpha1)), method='recursive', init=c(0,0))
    
    super.smoother.filter(HP)
}

# My Stochastic Indicator 2013 John F. Ehlers
roofing.stochastic.indicator  <- function(x, lookback = 20) {
    Filt = roofing.filter(x)
    
    HighestC = runMax(Filt, lookback)
        HighestC[1:lookback] = as.double(HighestC[lookback])
    LowestC = runMin(Filt, lookback)
        LowestC[1:lookback] = as.double(LowestC[lookback])
    Stoc = (Filt - LowestC) / (HighestC - LowestC)

    super.smoother.filter(Stoc)
}

First let’s plot the John’s stochastic oscillator vs traditional one.

###############################################################################
# Load Systematic Investor Toolbox (SIT)
# https://systematicinvestor.wordpress.com/systematic-investor-toolbox/
###############################################################################
setInternet2(TRUE)
con = gzcon(url('http://www.systematicportfolio.com/sit.gz', 'rb'))
    source(con)
close(con)

    #*****************************************************************
    # Load historical data
    #******************************************************************   
    load.packages('quantmod')  
    
    tickers = spl('DG')
    data = new.env()
    getSymbols(tickers, src = 'yahoo', from = '1970-01-01', env = data, auto.assign = T)   
        for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)
    bt.prep(data)

    #*****************************************************************
    # Setup
    #*****************************************************************
    prices = data$prices  
   
    models = list()
        
    # John Ehlers Stochastic
    stoch = roofing.stochastic.indicator(prices)
    
    # 14 Day Stochastic
    stoch14 = bt.apply(data, function(x) stoch(HLC(x),14)[,'slowD'])
		
    #*****************************************************************
    # Create plots
    #******************************************************************           
    dates = '2011:10::2012:9'
           
    layout(1:3)
    
    plota(prices[dates], type='l', plotX=F)
    plota.legend('DG')
    
    plota(stoch[dates], type='l', plotX=F)
        abline(h = 0.2, col='red')
        abline(h = 0.8, col='red')
    plota.legend('John Ehlers Stochastic')

    plota(stoch14[dates], type='l')
        abline(h = 0.2, col='red')
        abline(h = 0.8, col='red')
    plota.legend('Stochastic')

Next let’s code the trading rules as in Figure 6 and 8 in the John Ehlers paper: Predictive Indicators for Effective Trading Strategies

    #*****************************************************************
    # Code Strategies
    #*****************************************************************
    # Figure 6: Conventional Wisdom is to Buy When the Indicator Crosses Above 20% and 
    # To Sell Short when the Indicator Crosses below 80%
    data$weight[] = NA
        data$weight[] = iif(cross.up(stoch, 0.2), 1, iif(cross.dn(stoch, 0.8), -1, NA))
    models$post = bt.run.share(data, clean.signal=T, trade.summary=T)

    data$weight[] = NA
        data$weight[] = iif(cross.up(stoch, 0.2), 1, iif(cross.dn(stoch, 0.8), 0, NA))
    models$post.L = bt.run.share(data, clean.signal=T, trade.summary=T)
    
    data$weight[] = NA
        data$weight[] = iif(cross.up(stoch, 0.2), 0, iif(cross.dn(stoch, 0.8), -1, NA))
    models$post.S = bt.run.share(data, clean.signal=T, trade.summary=T)
        
    # Figure 8: Predictive Indicators Enable You to Buy When the Indicator Crosses Below 20% and 
    # To Sell Short when the Indicator Crosses Above 80%
    data$weight[] = NA
        data$weight[] = iif(cross.dn(stoch, 0.2), 1, iif(cross.up(stoch, 0.8), -1, NA))
    models$pre = bt.run.share(data, clean.signal=T, trade.summary=T)

    data$weight[] = NA
        data$weight[] = iif(cross.dn(stoch, 0.2), 1, iif(cross.up(stoch, 0.8), 0, NA))
    models$pre.L = bt.run.share(data, clean.signal=T, trade.summary=T)

    data$weight[] = NA
        data$weight[] = iif(cross.dn(stoch, 0.2), 0, iif(cross.up(stoch, 0.8), -1, NA))
    models$pre.S = bt.run.share(data, clean.signal=T, trade.summary=T)
    	
    #*****************************************************************
    # Create Report
    #****************************************************************** 		    
    strategy.performance.snapshoot(models, T)

Look’s like most profit comes from the long side.

Finally let’s visualize the timing of the signal for these strategies:

john.ehlers.custom.strategy.plot <- function(data, models, name, 
	main = name,
	dates = '::',
	layout = NULL		# flag to indicate if layout is already set	
) {
    # John Ehlers Stochastic
    stoch = roofing.stochastic.indicator(data$prices)
    	
    # highlight logic based on weight
    weight = models[[name]]$weight[dates]
    	col = iif(weight > 0, 'green', iif(weight < 0, 'red', 'white'))
    plota.control$col.x.highlight = col.add.alpha(col, 100)
    	highlight = T
   
    if(is.null(layout)) layout(1:2)
    	 
    plota(data$prices[dates], type='l', x.highlight = highlight, plotX = F, main=main)
    plota.legend('Long,Short,Not Invested','green,red,white')
    	
    plota(stoch[dates], type='l', x.highlight = highlight, plotX = F, ylim=c(0,1))        	
       	col = col.add.alpha('red', 100)
        abline(h = 0.2, col=col, lwd=3)
        abline(h = 0.8, col=col, lwd=3)
    plota.legend('John Ehlers Stochastic')        
}

		
    layout(1:4, heights=c(2,1,2,1))

    john.ehlers.custom.strategy.plot(data, models, 'post.L', dates = '2013::', layout=T,
        main = 'post.L: Buy When the Indicator Crosses Above 20% and Sell when the Indicator Crosses Below 80%')

    john.ehlers.custom.strategy.plot(data, models, 'pre.L', dates = '2013::', layout=T,
	main = 'pre.L: Buy When the Indicator Crosses Below 20% and Sell when the Indicator Crosses Above 80%')

As a homework, I encourage you to compare trading the John’s stochastic oscillator vs traditional one.

To view the complete source code for this example, please have a look at the john.ehlers.filter.test() function in bt.test.r at github.

Longer-history back-tests

One of the important steps of evaluating new trading idea or strategy is to see how it behaved historically (i.e. create back-test and examine the equity curve in different economic and market conditions)

However, creating a long back-test is usually problematic because most ETFs do not have a long price history. One way to alleviate the short price history is to find an appropriate proxy asset/index to extend ETF’s price history.

For example in my last post, Update: Extending Commodity time series, I showed howyou can extend the price history of commodity ETFs with Thomson Reuters / Jefferies CRB Index.

I created two functions in data.proxy.r at github to help visualize and compute statistics for potential proxy time series:

• proxy.test function will create a plot and compute summary statistics over the common period for all given series
• proxy.overlay.plot function will create an overlay plot for all given series (a plot where all shorter series are overlayed over the longest one)

Following is an example of using these functions for commodity ETFs and Thomson Reuters / Jefferies CRB Index.

###############################################################################
# Load Systematic Investor Toolbox (SIT)
# https://systematicinvestor.wordpress.com/systematic-investor-toolbox/
###############################################################################
setInternet2(TRUE)
con = gzcon(url('http://www.systematicportfolio.com/sit.gz', 'rb'))
    source(con)
close(con)

    #*****************************************************************
    # Load historical data
    #******************************************************************   
    load.packages('quantmod')  
    
    tickers = spl('GSG,DBC')
    data = new.env()
    getSymbols(tickers, src = 'yahoo', from = '1970-01-01', env = data, auto.assign = T)   
        for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)
      
    # "TRJ_CRB" file was downloaded from the http://www.jefferies.com/Commodities/2cc/389
    # for "TRJ/CRB Index-Total Return"
    temp = extract.table.from.webpage( join(readLines("TRJ_CRB")), 'EODValue' )
    temp = join( apply(temp, 1, join, ','), '\n' )
    data$CRB_1 = make.stock.xts( read.xts(temp, format='%m/%d/%y' ) )
     
    # "prfmdata.csv" file was downloaded from the http://www.crbequityindexes.com/indexdata-form.php
    # for "TR/J CRB Global Commodity Equity Index", "Total Return", "All Dates"
    data$CRB_2 = make.stock.xts( read.xts("prfmdata.csv", format='%m/%d/%Y' ) )
                 
    #*****************************************************************
    # Compare
    #******************************************************************    
    proxy.test(data)    

    proxy.overlay.plot(data)

The Thomson Reuters / Jefferies CRB Index (CRB_1) looks a better fit for commodity ETFs.

Yahoo Finance contains a big collection of indices with long histories. Just to give you an example, following is a list of some indices for major asset classes:

• Equity Market
• Vanguard Total Stock Mkt (VTSMX)
• Vanguard Total Intl Stock (VGTSX)
• Vanguard 500 Index (VFINX)
• Vanguard Emerging Mkts (VEIEX)
• Fixed Income Market
• Vanguard Short-Term Treasury (VFISX)
• Vanguard Long-Term Treasury (VUSTX)
• Vanguard Total Bond Market (VBMFX)
• Vanguard High-Yield Corporate (VWEHX)
• PIMCO Emerging Markets Bond (PEBIX)
• Vanguard Inflation-Protected (VIPSX)
• PIMCO Total Return (PTTRX)
• Vanguard REIT (VGSIX)

Now let’s look at an example of extending real estate time ETFs:

    #*****************************************************************
    # Load historical data
    #******************************************************************   
    load.packages('quantmod')  
    
    tickers = spl('IYR,VGSIX,RWO')
    data = new.env()
    getSymbols(tickers, src = 'yahoo', from = '1970-01-01', env = data, auto.assign = T)   
        for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)
                 
    #*****************************************************************
    # Compare
    #******************************************************************    
    proxy.test(data)    

    proxy.overlay.plot(data)

The Vanguard REIT (VGSIX) Index looks as a good proxy for the Dow Jones U.S. Real Estate Index (IYR), while Dow Jones Global Real Estate Index (RWO) does behave differently.

So the next time you are faced with a task of running a long history back-test, please have a look at extending your short history assets with proxies.

Update: Extending Commodity time series

I showed an example of Extending Commodity time series back in 2012. Since then, the web site that I used to get the Thomson Reuters/Jefferies CRB Index data is no longer working. But there are a few alternatives:

• Thomson Reuters / Jefferies CRB Index.
  To get data, first select “TRJ/CRB Index-Total Return”, next click “See Chart”. There will be “Download to Spreadsheet” link beneath the chart.
• TR/J CRB Global Commodity Equity Index
  To get data, first select “TR/J CRB Global Commodity Equity Index”, “Total Return”, check “All Dates”, next click “OK” to download the data.

There is a difference between these two series. I did not found a reason for the discrepancy, but truthfully, I did not look very hard. So I welcome any insights or feedback.

The Thomson Reuters / Jefferies CRB Index time series is better fit for the DBC and GSG etfs.

Now let’s plot the CRB index alternatives and DBC and GSG Commodity ETFs side by side to visual check the similarity.

###############################################################################
# Load Systematic Investor Toolbox (SIT)
# https://systematicinvestor.wordpress.com/systematic-investor-toolbox/
###############################################################################
setInternet2(TRUE)
con = gzcon(url('http://www.systematicportfolio.com/sit.gz', 'rb'))
    source(con)
close(con)
 
    #*****************************************************************
    # Load historical data
    #******************************************************************    
    load.packages('quantmod')		
    tickers = spl('GSG,DBC')
    data = new.env()
    getSymbols(tickers, src = 'yahoo', from = '1970-01-01', env = data, auto.assign = T)    
        for(i in ls(data)) data[[i]] = adjustOHLC(data[[i]], use.Adjusted=T)
     
    # "TRJ_CRB" file was downloaded from the http://www.jefferies.com/Commodities/2cc/389
    # for "TRJ/CRB Index-Total Return"
    temp = extract.table.from.webpage( join(readLines("TRJ_CRB")), 'EODValue' )
	temp = join( apply(temp, 1, join, ','), '\n' )
	data$CRB_1 = make.stock.xts( read.xts(temp, format='%m/%d/%y' ) )
	
    # "prfmdata.csv" file was downloaded from the http://www.crbequityindexes.com/indexdata-form.php
    # for "TR/J CRB Global Commodity Equity Index", "Total Return", "All Dates"
    data$CRB_2 = make.stock.xts( read.xts("prfmdata.csv", format='%m/%d/%Y' ) )
        	    
    bt.prep(data, align='remove.na')
	    
    #*****************************************************************
    # Compare
    #******************************************************************    	
    plota.matplot(scale.one(data$prices))