MATLAB Examples

# Algorithmic Trading with MATLAB®: Simple Lead/Lag EMA

This demo is an introduction to using MATLAB to develop and test a simple trading strategy using an exponential moving average.

## Load in some data (Excel)

Bund is a German bond future and data is sampled daily

```data = xlsread('BundDaily.xls'); testPts = floor(0.8*length(data(:,5))); BundClose = data(1:testPts,5); BundCloseV = data(testPts+1:end,5); ```

## Develop a simple lead/lag technical indicator

We'll use two exponentially weighted moving averages

```[lead,lag]=movavg(BundClose,5,20,'e'); plot([BundClose,lead,lag]), grid on legend('Close','Lead','Lag','Location','Best') ``` Develop a trading signal and performance measures. We'll assume 250 trading days per year.

```s = zeros(size(BundClose)); s(lead>lag) = 1; % Buy (long) s(lead<lag) = -1; % Sell (short) r = [0; s(1:end-1).*diff(BundClose)]; % Return sh = sqrt(250)*sharpe(r,0); % Annual Sharpe Ratio ```

Plot results

```ax(1) = subplot(2,1,1); plot([BundClose,lead,lag]); grid on legend('Close','Lead','Lag','Location','Best') title(['First Pass Results, Annual Sharpe Ratio = ',num2str(sh,3)]) ax(2) = subplot(2,1,2); plot([s,cumsum(r)]); grid on legend('Position','Cumulative Return','Location','Best') linkaxes(ax,'x') ``` ## Sidebar: Single moving average

The case of a single moving average. We can use this function to do a single moving average by setting first parameter to 1.

```annualScaling = sqrt(250); leadlag(BundClose,1,20,annualScaling) ``` ## Sidebar: Best parameter

Perform a parameter sweep to identify the best setting.

```sh = nan(100,1); for m = 2:100 [~,~,sh(m)] = leadlag(BundClose,1,m); end [~,mxInd] = max(sh); leadlag(BundClose,1,mxInd,annualScaling) ``` ## Estimate parameters over a range of values

Return to the two moving average case and identify the best one.

```sh = nan(100,100); tic for n = 1:100 for m = n:100 [~,~,sh(n,m)] = leadlag(BundClose,n,m,annualScaling); end end toc ```
```Elapsed time is 2.900633 seconds. ```

Plot results

```figure surfc(sh), shading interp, lighting phong view([80 35]), light('pos',[0.5, -0.9, 0.05]) colorbar ``` Plot best Sharpe Ratio

```[maxSH,row] = max(sh); % max by column [maxSH,col] = max(maxSH); % max by row and column leadlag(BundClose,row(col),col,annualScaling) ``` ## Evaluate performance on validation data

```leadlag(BundCloseV,row(col),col,annualScaling) ``` ## Include trading costs

We'll add the trading cost associated with the bid/ask spread. This will get us closer to the actual profit we could expect. As an exercise, you should extend this to account for additional trading costs and slippage considerations.

```cost=0.01; % bid/ask spread range = {1:1:120,1:1:120}; annualScaling = sqrt(250); llfun =@(x) leadlagFun(x,BundClose,annualScaling,cost); tic [maxSharpe,param,sh,vars] = parameterSweep(llfun,range); toc figure surfc(vars{1},vars{2},sh), shading interp, lighting phong title(['Max Sharpe Ratio ',num2str(maxSharpe,3),... ' for Lead ',num2str(param(1)),' and Lag ',num2str(param(2))]); view([80 35]), light('pos',[0.5, -0.9, 0.05]) colorbar figure leadlag(BundCloseV,row(col),col,annualScaling,cost) ```
```Elapsed time is 4.452227 seconds. ```  ## Determine best trading frequency (considering intraday)

Load in 1-minute data and break into test/validation data sets

```close all load bund1min testPts = floor(0.8*length(data)); BundClose = data(1:testPts,4); BundCloseV = data(testPts+1:end,4); cost=0.01; % bid/ask spread ```

Best Lead/Lag model for minute data with frequency consideration. Use parallel computing to speed up the computations (parfor in leadlagFun)

```type leadlagFun ```
```function sh = leadlagFun(x,data,scaling,cost) % define leadlag to accept vectorized inputs and return only sharpe ratio %% % Copyright 2010, The MathWorks, Inc. % All rights reserved. [row,col] = size(x); sh = zeros(row,1); t = length(data); x = round(x); if ~exist('scaling','var') scaling = 1; end if ~exist('cost','var') cost = 0; end % run simulation parfor i = 1:row if x(i,1) > x(i,2) sh(i) = NaN; %elseif x(i,1) > t || x(i,2) > t %sh(i) = NaN; else if col > 2 tindex = 1:x(i,3):t; % calculate scaling parameter for time sampling sc = sqrt(scaling^2 / x(i,3)); else tindex = 1:t; sc = scaling; end [~,~,sh(i)] = leadlag(data(tindex), x(i,1), x(i,2),sc,cost); end end ```

Use my the cores on my laptop (a quadcore with hyperthreading, so 8 virtual cores).

```matlabpool local 8 ```
```Destroying 1 pre-existing parallel job(s) created by matlabpool that were in the finished or failed state. Starting matlabpool using the 'local' configuration ... connected to 8 labs. ```

Perform the parameter sweep

```seq = [1:20 10:10:100]; ts = [1:4 5:5:55 60:10:180 240 480]; range = {seq,seq,ts}; annualScaling = sqrt(250*11*60); llfun =@(x) leadlagFun(x,BundClose,annualScaling,cost); tic [~,param,sh,xyz] = parameterSweep(llfun,range); toc leadlag(BundClose(1:param(3):end),param(1),param(2),... sqrt(annualScaling^2/param(3)),cost) xlabel(['Frequency (',num2str(param(3)),' minute intervals)']) ```
```Elapsed time is 74.028335 seconds. ``` Plot iso-surface

```figure redvals = 1.2:0.1:1.9; yelvals = 0.3:0.1:1; bluevals=0:0.1:0.4; isoplot(xyz{3},xyz{1},xyz{2},sh,redvals,yelvals,bluevals) set(gca,'view',[-21, 18],'dataaspectratio',[3 1 3]) grid on, box on % labels title('Iso-surface of Sharpes ratios.','fontweight','bold') zlabel('Slow Mov. Avg.','Fontweight','bold'); ylabel('Fast Mov. Avg.','Fontweight','bold'); xlabel('Frequency (minutes)','Fontweight','bold'); ``` Note that the lag of 100 is on the boundary of our parameter sweep, let's extend the search a bit more. I ran this earlier and the max is around 30 minutes, so we'll narrow our sweep (for time considerations).

```seq = [1:20 300:1:400]; ts = 25:50; range = {seq,seq,ts}; annualScaling = sqrt(250*11*60); llfun =@(x) leadlagFun(x,BundClose,annualScaling,cost); tic [maxSharpe,param,sh,xyz] = parameterSweep(llfun,range); toc param %#ok<NOPTS> leadlag(BundClose(1:param(3):end),param(1),param(2),... sqrt(annualScaling^2/param(3)),cost) xlabel(['Frequency (',num2str(param(3)),' minute intervals)']) ```
```Elapsed time is 100.845390 seconds. param = 10 394 29 ``` ## Best performer on validation data

This is the result if we applied it to the remaining 20% (validation set) of the data.

```leadlag(BundCloseV(1:param(3):end),param(1),param(2),... sqrt(annualScaling^2/param(3)),cost) xlabel(['Frequency (',num2str(param(3)),' minute intervals)']) ``` Let's now add an RSI indicator and see if we can do better (AlgoTradingDemo2.m).