Is there a quicker algorithm for calculating 'drifted' portfolio weights? (R, Pandas/NumPy, MATLAB)

Question

'Sup, QuantSX.

BLOT (Bottom Line On Top): Is there a nice clean algorithm for rapidly calculating portfolio weight drift? In R, Python or MATLAB - I'm not fussed which.

Details

I'm in the final stages of extending a stress testing script (in MATLAB) that I've used for a while, and I've hit a performance wall on what should be a pretty easy updating function.

The guts of the script is a medium-sized Monte Carlo - 1000 iterations of 240 x n[i] for ~70 portfolios, where n[i] is the number of portfolio components in the ith portfolio (n[i] ranges from ~5 to ~25). Gaussian copula, so generates the pseudo-data super fast.

In the past the output has just been presented under a monthly-rebalance assumption, so the pseudo-observations for the portfolio itself are just

$w\cdot R = \sum_i\bar{w}_i \times R_{i,t}$

where $\bar{w}_i$ are the original-and-target portfolio weights. And of course that's also super-quick - O(msec). In fact it's the quickest module of the whole stress-testing framework, which completes in ~35 seconds per portfolio.

Part of the extension is to calculate the portfolio path under a range of alternative rebalancing assumptions - full-drift (i.e., no rebalancing, buy-and-hold); quarterly; annual; and 'threshold' by asset class.

In testing the full-drift version, it's taking 2 seconds per iteration to calculate the drifted weights; that seems a very long time for 240 rows of data.

And 70 lots of 1000 × ~2 seconds is a long time.

The weight transition equation is really straightforward -

$w_{i,t+1} = w_{i,t}\cdot\left(\frac{1+R_{i,t}}{1+R_{p,t}}\right)$

$w_{i,0} = \bar{w}_i$

The drifted weight depends on the portfolio return, which is the weighted sum of the component returns in period t -

$R_{p,t} = \sum_i(w_{i,t}\cdot R_{i,t}) = w'̌\cdot R$

So far, no rocket surgery.

To do this in MATLAB, the only way I could think of was to initialise an aggregator, and to accumulate the recalculated weights (and as a byproduct, the 'drifted' portfolio returns).

So my function is presented below (although my function's properly indented)

function[outPVals, outWts, wtDiff] = DriftWts(inTable, inWts, outName) 
%% Calculates Portfolio Value and Component Weights over time if weights drift.
xT = inTable;
xP = inTable;
xT.Properties.DimensionNames{2} = 'Weights';
xT.Weights = zeros(size(xT)); %% pre-allocate weights table
xD = xT;                      %% pre-allocate weight-delta table 
tBeg = min(xT.MthYr);
xT(tBeg,:).Weights = transpose(inWts); %% initialises weights at current AA
xR = inTable(tBeg,:).Returns/100;
xP.(outName)(tBeg,:) = 100*(xR * inWts); 
xT = sortrows(xT,'MthYr','ascend');
for mmm = 2:length(xT.MthYr)
    m = xT.MthYr(mmm);
    m_l = eomdate(m - calmonths(1));
    w_l = xT(m_l,:).Weights;
    R_l = inTable(m_l,:).Returns/100;
    xR = inTable(m,:).Returns/100;
    Rp_l = w_l * transpose(R_l);
    w_0 = w_l .* (1 + R_l)/( 1 + Rp_l);
    xT(m,:).Weights = w_0;
    xD(m,:).Weights = w_0 - w_l;
    xP.(outName)(m,:) = 100*(w_0 * transpose(xR));
 end
 outWts = xT;
 outPVals = xP(:,outName);
 wtDiff = xD;
 end

I considered an alternative - i.e., generating cumprod(1 +$R_{i,t}$ /100) of the component returns, weighting them by $\bar{w}_i$ at t=0, summing them, then dividing everything through by the row sum.

It's not obvious to me that this alternative would be quicker (but I have no idea what goes on under the hood in cumprod(), so I'll try it today)

Answer 1

OK, so the 'alternative' mechanism (cumprod, weight, divide by row sum) is O(msec) and produces identical results.

• Old stupid method: 1.705sec for a 200×9 table.
• New sexy method: 0.004sec for the same table.

And it's only 4 lines of code.

T01 = cumprod(1+inTable/100,'reverse');
T02 = inWts' .* T01;
R01 = sum(T02, 2);
W01 = T02 ./ R01;