FindMinimum and least square minimization
- To: mathgroup at smc.vnet.net
- Subject: [mg30703] FindMinimum and least square minimization
- From: "Dr J. Satherley" <js1 at liverpool.ac.uk>
- Date: Sat, 8 Sep 2001 02:23:18 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
Dear All
I have a complicated non-linear least square minimization problem. It is to
fit numerical data to a function that has to be solved using NDSolve at each
iteration. To do this I have written a function to compute the sum of
squares between each data point and the value of the function at this point.
I then use FindMinimum to find the values of the parameters which minimise
the sum of squares. Mathematica does this very slowly and I want to find a
more efficient way to do the calculation. To that end I have worked on a
simple example to assist me with finding improvements, the main one of which
is to supply the partial derivatives of the function with respect to each
parameter. However, the example leaves me a little perplexed and I wonder if
anyone out there can enlighten me on the points I raise below.
Here's the example:
First the inital equation:
y[x_] := 2.31 x^3 + 1.2 x^2 + 3.42 Sin[x] + 131.56
Generate some data using this equation:
xydata = Table[{u, y[u]}, {u, -10, 10, 0.1}]
Write the function to calculate the sum of squares. I'm simply trying to
find the parameters in the above equation:
Q[a_?NumberQ,b_?NumberQ,c_?NumberQ,d_?NumberQ]:=Module[{output,fn,temp},
fn=#[[2]]-(a #[[1]]^3+b #[[1]]^2+c Sin[#[[1]]]+d )&;
temp=fn/@xydata;
output=temp.temp;
total=totala+totalb+totalc+totald;
Print[output," ",++n," ",a," ",b," ",c," ",d," ",totala," ",totalb,"
",totalc," ",totald," ",total];(*this prints out what is happening at each
iteration*)
output
];
Define the partial derivatives of Q wrt each parameter:
dQda[a_, b_, c_, d_] := totala = Plus @@ (2*(#1[[2]] - (a*#1[[1]]^3 +
b*#1[[1]]^2 + c*Sin[#1[[1]]]+d))*#1[[1]]^3 & ) /@ xydata;
dQdb[a_, b_, c_, d_] := totalb = Plus @@ (2*(#1[[2]] - (a*#1[[1]]^3 +
b*#1[[1]]^2 + c*Sin[#1[[1]]]+d))*#1[[1]]^2 & ) /@ xydata;
dQdc[a_, b_, c_, d_] := totalc = Plus @@ (2*(#1[[2]] - (a*#1[[1]]^3 +
b*#1[[1]]^2 + c*Sin[#1[[1]]]+d))*Sin[#1[[1]]] & ) /@ xydata;
dQdd[a_, b_, c_, d_] := totald = Plus @@ (2*(#1[[2]] - (a*#1[[1]]^3 +
b*#1[[1]]^2 + c*Sin[#1[[1]]]+d)) & ) /@ xydata;
Derivative[1, 0, 0, 0][Q] := dQda;
Derivative[0, 1, 0, 0][Q] := dQdb;
Derivative[0, 0, 1, 0][Q] := dQdc;
Derivative[0, 0, 0, 1][Q] := dQdd;
Run FindMinimum:
MemoryInUse[]
n=0;result2=FindMinimum[Q[a,b,c,d],{a,2},{b,2},{c,2},{d,2},MaxIterations->50
]
MemoryInUse[]
These are the points I've noted when running these functions:
1. I was expecting the convergence to be rather rapid compared to giving 2
starting values to FindMinimum. However, it is only marginally quicker -
maybe 150 iterations instead of 220. Is this to be expected? Or have I not
formulated the problem correctly?
2. I was expecting the sum of the 4 parital derivative functions to approach
zero at the convergence point. However, it was not as close as I would have
thought - for example only around 0.007.
3. Reaching the correct solution is more sensitive to the choice of starting
values when using FindMinimum together with the partial derivatives. Using
FindMinimum with 2 starting values does take more iterations but reaches a
solution with a smaller sum of squares (4 orders of magnitude less).)
4. I have noticed that Mathematica uses a vast amount of memory when I've
performed my actually problem. It uses up all the available RAM (256MB on my
system) and them runs off the harddisk using the swap file. That is why I've
included the MemoryInUse before and after running FindMinimum to monitor the
memory use. Even with my simple example memory is not returned for use. I'm
using Mathematica 4.0.1 on Windows98. Is there anything I can do to fix
this problem? On my real problem after running the function a couple of
times it slows dramatically once the harddisk is accessed.
I'd be grateful for any comments or suggestions related to these
observations particularly those that may reduce the number of iterations and
the problem about the memory use.
Cheers
John Satherley