| Author |
Comment/Response |
Jeremy Jones
|
12/09/11 08:48am
I use LinearModelFit to get a polynomial fit to some data (I choose the degree of the polynomial based on AIC). I then look for turning points using NSolve. But when I substitute the answers back into the derivative, some of them are not very close to zero. Here is my code.
data = {43.46, 43.63, 43.6, 42.89, 42.85, 42.62, 43.44, 43.58, 43.43, 43.06, 43.52, 43., 42.03, 41.42, 42.3, 41.57, 40.15, 40.06, 38.12, 36.6, 37.05, 38.4, 37.93, 38.21, 39.85, 39.43, 39.82, 39.1, 37.5, 37.53, 38.21, 38.21, 38.61, 38.64, 39.15, 39.35, 39.74, 39.87, 39.04, 37.7, 36.88, 38.23, 38.17, 37.91, 36.45, 37.29, 36.9, 37.63, 38.23, 37.6, 36.8, 37.13, 35.63, 34.55, 33.95, 35.35, 35.52, 35.04, 35.02, 34.15, 33.86, 35.13, 36.3, 37.2, 37.1, 37.39, 37.08, 37.64, 36.86, 37.65, 36.4, 36.4, 35.48, 35.7, 36.85, 37.17, 36.94, 38.35, 38.69, 37.8, 36.77, 36.55, 36.54, 37.95, 37.7, 37.75, 38.33, 37.48, 37.8, 37.5, 36.88, 36.65, 37.04, 36.13, 35.9, 35.6, 34.51, 34.53, 34.05, 34.85};
xs = Range[Length[data]];
ys = data;
pairs = MapThread[{#1, #2} &, {xs, ys}];
order = 24;
vars = x^(# - 1) & /@ Range[order + 1];
lm = LinearModelFit[pairs, vars, x];
fn = Function[x, Evaluate[lm["BestFit"]]];
xmin = First[xs] - 0.5;
xmax = Last[xs] + 0.5;
tps = NSolve[fn'[x] == 0 && xmin < x < xmax, x];
Substituting the answers back into fn'[x] gives
fn'[x] /. tps
{6.55301*10^-15, -8.32358*10^-14, -2.49786*10^-13, 1.16063*10^-12, -4.33076*10^-12, 1.16415*10^-10, -9.31323*10^-10, -7.45058*10^-9, 1.49012*10^-8, 7.62939*10^-6, -7.62939*10^-6, 0.000167847, 0.000366211, 0.0012207, -0.00146484, 0.0200195, 0.0200195}
The values get rather large towards the end. Also, I check the nature of the turning points with the second derivative, and I get
fn''[x] /. tps
{-4.16667, 0.580779, -0.164657, 0.106349, -0.165436, 0.196445, -0.0844069, 0.0672132, -0.0899606, 0.212267, -0.199458, 0.18417, -0.172743, 0.180859, -0.292915, -0.00140381, -0.00140381}
The last few are all negative, but they should alternate between positive and negative.
I'm not sure what the problem is, it might have something to do with small coefficients in the polynomial. Help would be appreciated. Thanks.
URL: , |
|