Re: Best linear Fit to slope data with Fixed starting point/value.
- To: mathgroup at smc.vnet.net
- Subject: [mg63592] Re: Best linear Fit to slope data with Fixed starting point/value.
- From: "Ray Koopman" <koopman at sfu.ca>
- Date: Fri, 6 Jan 2006 05:24:40 -0500 (EST)
- References: <dpg13b$pmj$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Lea Rebanks wrote:
> Given the list below, (quite large - 4516 points), is when plotted an
> almost perfect straight line with noise.
>
> I want to plot the best linear fit (i.e. m*x + b) to this data, however
> I must retain the original first value (151.0993767999595).
>
> Obviously I could do this manually, but was wondering if there was a
> more accurate & efficient way using Mathematica.
>
> Many thanks for your attention.
>
> Best Regards - Lea Rebanks...
>
>
> Data below:-
>
> { (* data omitted *) }
Why are there only 1929 points, instead of 4516 ?
In[1]:= n = Length[y = { (* data omitted *) }]
Out[1]= 1929
This data set requires careful analysis. The purpose of the many
"InputForm"s in the following is to force all the digits to print.
f is the usual unconstrained least-squares fit.
r is the corresponding residual vector,
followed by some descriptive statistics.
Mean, Min, & Max ought to be zero, negative, & positive, respectively.
In[2]:= {f = Fit[y,{1,x},x], r = y - Table[f,{x,n}];
r[[1]], Mean@r, Min@r, Max@r, r.r} //InputForm
Out[2]//InputForm=
{151.0405970247262 + 0.05877977523343556*x,
-1.4210854715202004*^-13, -1.6650789957127822*^-13,
-2.2737367544323206*^-13, -1.4210854715202004*^-13,
5.425222374968153*^-23}
f0 is a least-squares fit, constrained to pass thru {1,y[[1]]}.
r0 is the corresponding residual vector,
followed by the same descriptive statistics as for r.
r0[[1]] ought to be zero, and r0.r0 ought to be greater than r.r .
In[3]:= {f0 = y[[1]] + Fit[y-y[[1]],{x-1},x], r0 = y - Table[f0,{x,n}];
r0[[1]], Mean@r0, Min@r0, Max@r0, r0.r0} //InputForm
Out[3]//InputForm=
{151.0993767999595 + 0.05877977523343556*(-1 + x),
0., -2.1305231641453703*^-14,
-1.1368683772161603*^-13, 0.,
1.5832753912147795*^-24}
f & f0 differ only in their intercepts.
In[4]:= Simplify[f - f0] //InputForm
Out[4]//InputForm= 1.4519635493925875*^-13 + 0.*x
Doing everything in double precision gives better results.
In[5]:= yy = SetPrecision[y,2*$MachinePrecision];
{ff = Fit[yy,{1,x},x], rr = yy - Table[ff,{x,n}];
rr[[1]], Mean@rr, Min@rr, Max@rr, rr.rr} //N//InputForm
Out[6]//InputForm=
{151.04059702472605 + 0.05877977523343554*x,
2.0177667765604405*^-15, 0.,
-4.916104707204294*^-14, 3.285806817257583*^-14,
2.475826245403202*^-25}
In[7]:= {ff0 = yy[[1]] + Fit[yy-yy[[1]],{x-1},x],
rr0 = yy - Table[ff0,{x,n}];
rr0[[1]], Mean@rr0, Min@rr0, Max@rr0, rr0.rr0} //N//InputForm
Out[7]//InputForm=
{151.0993767999595 + 0.05877977523343554*(-1. + x),
0., -5.048340522117773*^-16,
-4.833186368058629*^-14, 3.3630752001712404*^-14,
2.4954757604292457*^-25}
In[8]:= Simplify[ff - ff0] //N//InputForm
Out[8]//InputForm= -2.019336208847109*^-15 + 1.5694322866687375*^-18*x
y appears to be quantized,
with more noise at the top than at the bottom.
In[9]:= y1 = ListConvolve[{1,-1},y];
ListPlot[y1,PlotRange->All,Frame->True,Axes->None];
g = {#[[1]],Length@#}& /@ Split@Reverse@Sort@y1;
TableForm[{InputForm@#[[1]],#[[2]]}&/@g]
Out[12]//TableForm= 0.058779775233460896 280
0.058779775233432474 1577
0.05877977523340405 70
0.05877977523337563 1
In[13]:= y2 = ListConvolve[{1,-1},y1];
ListPlot[y2,PlotRange->All,Frame->True,Axes->None];
h = {#[[1]],Length@#}& /@ Split@Reverse@Sort@y2;
TableForm[{InputForm@#[[1]],#[[2]]}&/@h]
Out[16]//TableForm= 8.526512829121202*^-14 1
5.684341886080802*^-14 70
2.842170943040401*^-14 192
0. 1402
-2.842170943040401*^-14 191
-5.684341886080802*^-14 70
-8.526512829121202*^-14 1
In[17]:= ( g[[All,1]] - g[[2,1]] ) / h[[3,1]] //InputForm
Out[17]//InputForm= {1., 0., -1., -2.}
In[18]:= {y===Sort@y, m = (n+1)/2, y[[m]], Mean@y - y[[m]]} //InputForm
Out[18]//InputForm= {True, 965, 207.76308012499135, 0.}
In[19]:= z = y - y[[m]]; z[[m]]
k = {#[[1]],Length@#}& /@ Split@Reverse@Sort@Drop[z + Reverse@z, m];
TableForm[{InputForm@#[[1]],#[[2]]}&/@k]
Out[19]= 0.
Out[21]//TableForm= 2.842170943040401*^-14 95
0. 766
-2.842170943040401*^-14 76
-5.684341886080802*^-14 27
In[22]:= k[[1,1]] === h[[3,1]]
Out[22]= True
Is 2.842170943040401*^-14 some sort of quantum number?