Help: Using fminsearch to calculate fractional derivatives

Hi everyone,

I have been trying to incorporate caputo derivatives into my equation to calculate the optimal fractional derivatives and coefficients for the corresponding variables. I am using MATLAB'S inbuilt function fminsearch by minisizing the MSE to get the answer. However, my code isnt running when I tried to implement them on MATLAB. Can I get some enlightenment on how should I edit these code? (I have the caputo functions running perfectly with 5 parameters dy=caputo(y,t,alfa,y0,L), where

% y - the samples or function handle of the original function

% t - the time vector

% alpha - the fractional order

% y0 - the initial vector of signal and its integer-order derivatives

% L0 - the interpolation length

% dy - the Caputo derivative)

% Original Equation: Y=c(1)*X1+c(2)*X2+c(3)*X3+c(4)*x4+c(5)*X5+c(6)*X6+c(7)*X7+c(8)*X8+c(9)*X9
% Variable Matrix:
Y=GDP*10^11;
Y = Y(:);
X = [caputo(X1,t,c(10),91500, 1:53) caputo(X2/100*X1,t,c(11),2507100, 1:53) caputo(X3,t,c(12),8857716, 1:53) caputo(X4,t,c(13),1.4, 1:53) caputo(X5*10^9,t,c(14),7.49*10^9, 1:53) caputo(X6*10^9,t,c(15),5.18*10^9, 1:53) caputo(X7*10^9,t,c(16),3.32*10^9, 1:53) caputo(X8*10^10,t,c(17),1.88*10^10, 1:53) caputo(X5*10^9,t,c(18),7.49*10^9, 1:53)];
X = [caputo(X1,t,c(10),91500, 1:53)(:) caputo(X2/100*X1,t,c(11),2507100, 1:53)(:) caputo(X3,t,c(12),8857716, 1:53)(:) caputo(X4,t,c(13),1.4, 1:53)(:) caputo(X5*10^9,t,c(14),7.49*10^9, 1:53)(:) caputo(X6*10^9,t,c(15),5.18*10^9, 1:53)(:) caputo(X7*10^9,t,c(16),3.32*10^9, 1:53)(:) caputo(X8*10^10,t,c(17),1.88*10^10, 1:53)(:) caputo(X5*10^9,t,c(18),7.49*10^9, 1:53)(:)];
% Parameters: c1=c(1) c2=c(2) c3=c(3) cn=c(n) where n=[1,18]
Yfcn=@(c,X) X(:,1).*c(1)+X(:,2).*c(2)+X(:,3).*c(3)+X(:,4).*c(4)+X(:,5).*c(5)+X(:,6).*c(6)+X(:,7).*c(7)+X(:,8).*c(8)+X(:,9).*c(9)
SSECF=@(c) sum((Y - Yfcn(c,X)).^2)/53;
c0=-100*1e10+(2*100*1e10).*rand(1,18);             % Initial Parameter Estimates (-100*10^10,100*10^10)
options = optimset('MaxFunEvals', 1E7, 'MaxIter', 1E7);
[c,SSE] = fminsearch(SSECF,c0,options);

[J. Alex Lee]

not sure what you mean by "not running", but one thing is to consider using nsqnonlin or something similar; fminsearch will be ineffcieint for a 18-dimensional search

[Matt J]

@Yew Seng Tan There is a c in the lines that generates X

X = [caputo(X1,t,c(10),91500, 1:53) ...

This is not part of your objective function, even though you later use a different variable of the same name in your objective function SSECF. Where does this original c come from?

[Yew Seng Tan]

I'm trying to set up a fractional order derivatives for each of the variables X1-X9. c(10) to c(18) represent the order of fractional derivatives for each of the variables. I saw a paper stating that this alphas could be generated using the fminsearch function via MATLAB.

Answer 1

[Matt J]

Your optimization problem appears to be an unconstrained, linear least squares minimization. I don't see why you need an iterative solver at all. You could just solve with,

c=X\Y;

[J. Alex Lee]

@Matt J, at first i thought so too, but c10~c18 appear also as inputs into the caputo() function, which presumably makes everything nonlinear.

[Matt J]

@J. Alex Lee If that's what's intended, it isn't what's implemented in the posted code. fminsearch is seeing X as constant numeric data, not a further function of c.

[J. Alex Lee]

good point, only Yew Seng can clarify intent.

[Yew Seng Tan]

Yes c10-c18 are inputs to caputo. How should I change the code to allow fminsearch to obtain a good value for the fractional order derivatives?

[Matt J]

First, you must create a function that returns X as a function of c.

Xfcn=@(c,otherstuff)....

Then you must write a model for Y that utilizes Xfcn.

Yfcn=@(c,otherstuff) Xfcn(c,otherstuff)*reshape(c(1:9),[],1);

Then, you should give this model function to lsqcurvefit() with xdata=otherstuff and ydata=Y.

.

[Yew Seng Tan]

Thanks I will try this out. Can I get the coefficents from lsqcurvefit() too?

[Matt J]

It will give you c.

[Yew Seng Tan]

Thanks for the reply! Appreciate it! I will try this out now!