help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

leasqr problem in the partial derivative


From: Doug Stewart
Subject: leasqr problem in the partial derivative
Date: Thu, 21 Sep 2017 20:08:05 -0400

I went to the examples of leasqr and copied example #2, then adjusted it for my data.

% Define functions
 leasqrfunc = @(x, p) p(1) * exp (-p(2) * x ) +p(3);
 leasqrdfdp = @(x, f, p, dp, func) [exp(-p(2)*x), -p(1)*x.*exp(-p(2)*x),0*x+1];

what I did was added +p(3) to the original function. I then calculated the partial diff with respect to p(3). this is " 1 ".
when I put a 1 in the 3rd position in leasqrdfdp, then the program does not work. but when I make the 3rd term to by 0*x + 1  then the program runs and gives the correct answer.

should just +1 work???

here is the program that works 


% Define functions
 leasqrfunc = @(x, p) p(1) * exp (-p(2) * x ) +p(3);
 leasqrdfdp = @(x, f, p, dp, func) [exp(-p(2)*x), -p(1)*x.*exp(-p(2)*x),0*x+1];

 % generate test data
 t = [0 2 4 6 13]';
 #p = [1; 0.1];
 #data = "" (t, p);
 data=""  1457   1426   1408   1387]';
 #rnd = [0.352509; -0.040607; -1.867061; -1.561283; 1.473191; ...
 #       0.580767;  0.841805;  1.632203; -0.179254; 0.345208];

 % add noise
 % wt1 = 1 /sqrt of variances of data
 % 1 / wt1 = sqrt of var = standard deviation
 wt1 = (1 + 0.01 * t) ./ sqrt (data); 
 #data = "" + 0.05 * rnd ./ wt1; 
#wt1=[1 1 1 1 1]' *1;
 % Note by Thomas Walter :
 %
 % Using a step size of 1 to calculate the derivative is WRONG !!!!
 % See numerical mathbooks why.
 % A derivative calculated from central differences need: s 
 %     step = 0.001...1.0e-8
 % And onesided derivative needs:
 %     step = 1.0e-5...1.0e-8 and may be still wrong

 F = leasqrfunc;
 dFdp = leasqrdfdp; % exact derivative
 % dFdp = dfdp;     % estimated derivative
 dp = [01; 0.1;1];
 pin = [160; 1.0; 1150]; 
 stol=0.000001; niter=50;
 minstep = [0.01; 0.00001;.1];
 maxstep = [010; 0.2; 1000];
 options = [minstep, maxstep];

 global verbose;
 verbose = 1;
 [f1, p1, kvg1, iter1, corp1, covp1, covr1, stdresid1, Z1, r21] = ...
    leasqr (t, data, pin, F, stol, niter, wt1, [], dFdp, options);
    
    f1
    p1
    xx=0:.1:30;
 yy=F(xx,p1);
 hold on
 plot(xx,yy)
    hold off
    




--
DASCertificate for 206392


reply via email to

[Prev in Thread] Current Thread [Next in Thread]