Help: Run function for subsections of data

Hey everyone,

i want to loop or index my data in such a way, that MatLab repeats my lsqcurvefit function for all data, which falls on the same day.

I have options data for 818912 calls/puts over a span of 6151 days.

E = grp2idx(SPX.date);
t = unique (E);
Nt = histc(E(:),t);

gives me a vector of 6151 values, counting how many options there were on each trading day. (Nt = [10;11;13;13;9.....])

In general, i want to determine 10 parameters, which minimize the estimation error of my model. My function looks like this:

fun = @(x,xdata) (x(1)+x(2).*xdata(:,2)) + (x(3)+x(4).*xdata(:,2)).*((x(5)+x(6).*xdata(:,2)).*...
    (xdata(:,1) -(x(7)+x(8).*xdata(:,2))) + sqrt((xdata(:,1) - (x(7)+x(8).*xdata(:,2)).^2 ...
    + (x(9)+x(10).*xdata(:,2)).^2)));

where

A = [log_moneyness maturity]; xdata = A; ydata = implied_volatility;

Implied_volatility, log_moneyness and maturity are all observed data from the market for each day.

To give an example, what they look like (each 818912x1)

implied_volatility = [0.0525854412050268;0.0605709606968685;0.0736100506516179;0.0534246390268888;0.0370060201846111;...]
    maturity = [44;44;44;44;44;72;.....]
    log_moneyness=[-0.404048810868278;-0.233367535419970;-0.0640683097545084;0.103871068092051;0.270472269348038;...]

I continued with:

    %Fit the model using the starting point from former parameter estimation
    %via SVI 
    x0 = [0.000202758890760171 0.000202758890760171 -0.249601998773279 -0.249601998773279 ...
    0.5 0.5 0.0114 0.0114 1.73380716030294e-05 1.73380716030294e-05];
    
    opt = optimoptions (@lsqcurvefit, 'StepTolerance', 1e-10);
    
    x = lsqcurvefit (fun,x0,xdata,ydata,[],[],opt);

This gives me 10 parameters for all of my data. What i now want to get is 10 parameters for each of my 6151 observation days, so a matrix of 10x6151.

What i got until now is:

    for j=1:6151
    for    i=Nt(j);
        fun = @(x,xdata) (x(1)+x(2).*xdata(j:j+i-1,2)) + (x(3)+x(4).*xdata(j:j+i-1,2)).*((x(5)+x(6).*xdata(j:j+i-1,2)).*...
            (xdata(j:j+i-1,1) -(x(7)+x(8).*xdata(j:j+i-1,2))) + sqrt((xdata(j:j+i-1,1) - (x(7)+x(8).*xdata(j:j+i-1,2)).^2 + (x(9)+x(10).*xdata(j:j+i-1,2)).^2)));
        x0 = [0.000202758890760171 0.000202758890760171 -0.249601998773279 -0.249601998773279 ...
            0.5 0.5 0.0114 0.0114 1.73380716030294e-05 1.73380716030294e-05];
        opt = optimoptions (@lsqcurvefit, 'StepTolerance', 1e-10);
        x = lsqcurvefit (fun,x0,xdata(j:i-1,:),ydata(j:i-1,:),[],[],opt);     
    end
    end

Unfortunately, this always end in : Index in position 1 exceeds array bounds (must not exceed 9).

Can anyone tell me or help me to receive the parameters for the first day of data, then the second and so on? Let me know if there is anything unclear.

All the Best

Kai

[Kai Koslowsky]

Update:

for i=1:6151
    kk = N(i)
    ii = N (i+1);
    xdata= xdata('kk':'ii',:);
    ydata= ydata('kk':'ii',:);
    
    fun = @(x,xdata) (x(1)+x(2).*xdata(:,2)) + (x(3)+x(4).*xdata(:,2)).*((x(5)+x(6).*xdata(:,2)).*...
        (xdata(:,1) -(x(7)+x(8).*xdata(:,2))) + sqrt((xdata(:,1) - (x(7)+x(8).*xdata(:,2)).^2 ...
        + (x(9)+x(10).*xdata(:,2)).^2)));
    
    %Fit the model using the starting point from former parameter estimation
    %via SVI 
    x0 = zeros(1,10);
    
    %opt = optimoptions (@lsqcurvefit, 'StepTolerance', 1e-10);
    opt = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');
    
    x = lsqcurvefit (fun,x0,xdata,ydata,[],[],opt);
    parameters(i,:)= x;       
end

This is what i have until now, but i only gives me zeros for my parameter matrix.

Answer 1

[Walter Roberson]

'kk':'ii'

Warning: Colon operands must be real scalars.

ans =

  1×0 empty char array

You are passing empty xdata down...

fun = @(x,xdata) (x(1)+x(2).*xdata(:,2)) + (x(3)+x(4).*xdata(:,2)).*((x(5)+x(6).*xdata(:,2)).*...
        (xdata(:,1) -(x(7)+x(8).*xdata(:,2))) + sqrt((xdata(:,1) - (x(7)+x(8).*xdata(:,2)).^2 ...
        + (x(9)+x(10).*xdata(:,2)).^2)));

%Fit the model using the starting point from former parameter estimation
%via SVI 
x0 = [0.000202758890760171 0.000202758890760171 -0.249601998773279 -0.249601998773279 ...
0.5 0.5 0.0114 0.0114 1.73380716030294e-05 1.73380716030294e-05];
    
%opt = optimoptions (@lsqcurvefit, 'StepTolerance', 1e-10);
opt = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');

for i=1:6151
    kk = N(i)
    ii = N (i+1);
    Xdata = xdata(kk:ii,:);
    Ydata = ydata(kk:ii,:);
    x = lsqcurvefit (fun, x0, Xdata, Ydata, [], [], opt);
    parameters(i,:) = x;
end

I am not sure what N is in this context, though.

[Walter Roberson]

Also I suspect you should be using

ii = N(i+1) - 1;

Just make sure that the last N value is 1 greater than the last row number.

[Kai Koslowsky]

Hello Mr. Roberson,

i used:

 N = cumsum(Nt);

in order to receive only the values of xdata and ydata, which are specific to one day. So that MatLab uses only the data for one day and tries to fit those into lsqcurvefit.

I do not quite understand why I need this:

 Xdata = xdata(kk:ii,:);
 Ydata = ydata(kk:ii,:);

Could you elaborate? I am just running your code now, and will let you know, when it is done. Thank you very much already!

[Walter Roberson]

You were using

    xdata= xdata('kk':'ii',:);

That tried to take the character vector ['k' 'k'] and take the range to the character vector ['i' 'i'] .

When you use the colon operator with a non-scalar value, MATLAB takes the first element. So 'kk':'ii' is equivalent to 'k':'i' . Which is the same as char(double('k'):double('i'))) which is char(107:105) . where 107 is the character code for 'k' and 105 is the character code for 'i' . But 107 is after 105 so the result is empty.

Thus, your code was equivalent to

    xdata= xdata([], :)

which reads out emptiness on the right hand side, and assigns it over the complete variable xdata . Afterwards, xdata is empty . You took your non-empty xdata, extracted nothing from it, and overwrote the entire data with that emptiness.

My replacement code

    Xdata = xdata(kk:ii,:);

uses the values of the variables named kk and ii to extract a range of rows from xdata, and writes it to a different variable, named Xdata (capital-X), leaving xdata (lower-case X) unchanged. Then later in the code, it passes Xdata (the extracted rows) into the curvefitting function. The next round, the original xdata (lower-case X) is still all there to have different rows extracted.

[Kai Koslowsky]

I added a row above so the cumsum function takes the correct values.

Thank you very much, you made my day/night.

Greetings and all the best from Germany.

[Kai Koslowsky]

Thanks for the great explanation! As you might have guessed, i am fairly new to MatLab.

[Kai Koslowsky]

Is there an easy way to use particleswarm instead of lsqcurvefit for this problem? I was just told to use it instead and have not found a solution yet. All the best.

[Walter Roberson]

You can construct a new layer of function which is

residue = @(x) sum(fun(x)-Ydata).^2)

and use an optimizer on that. This value cannot be negative for real-valued function results, and smallest value corresponds to best fit.

[Kai Koslowsky]

Hello Mr Roberson,

thanks for your answer. Great to have people like you! Here is what i did so far:

A = [log_moneyness maturity]; xdata = A; ydata = implied_volatility; 
fun = @(x,xdata) (x(1)+x(2).*xdata(:,2)) + (x(3)+x(4).*xdata(:,2)).*((x(5)+x(6).*xdata(:,2)).*...
        (xdata(:,1) -(x(7)+x(8).*xdata(:,2))) + sqrt((xdata(:,1) - (x(7)+x(8).*xdata(:,2)).^2 ...
        + (x(9)+x(10).*xdata(:,2)).^2))); %original function

parameters_forall = table2array(T_param (:,2:11)); %1415596x10 matrix

%implied_volatility is a 1415596x1 matrix with original values from market
%now i want to minimize the error of my estimate by:   
tDays = size(parameters_forall(:,1));
%%
    for i = 1:tDays
        x = parameters_forall(i,:);
        residue = @(x) (sum(fun(x,xdata)-implied_volatility(i)).^2);
        y = particleswarm(residue, 10);
        swarm_parameters(i,:) = y;  
    end

So i am currently running this code. It already gave me a y=1415596x1 matrix, and is still going. Is particleswarm optimizing each x variable at a time? For the first, it gave me:

    Optimization ended: relative change in the objective value 
    over the last OPTIONS.MaxStallIterations iterations is less than OPTIONS.FunctionTolerance.

Does the syntax make sense to you for what i am trying to achive? (adjust all 10 parameters for each data on each day). All the best.

[Kai Koslowsky]

Okay, this takes for ever to complete. I stopped it after 30 minutes and it only had 10x12 parameters. Is there any way to speed up this code? Or should i put this as a new question?

[Walter Roberson]

You already created a new question ;-)