Re: fixed points piecewise-linear fitting

[Juan Pablo Carbajal]

On Sat, Mar 17, 2012 at 3:19 PM, Sergei Steshenko <address@hidden> wrote:
>
>
>
>
> ----- Original Message -----
>> From: Juan Pablo Carbajal <address@hidden>
>> To: Sergei Steshenko <address@hidden>
>> Cc: "address@hidden" <address@hidden>
>> Sent: Saturday, March 17, 2012 2:17 PM
>> Subject: Re: fixed points piecewise-linear fitting
>>
>> On Sat, Mar 17, 2012 at 12:43 PM, Sergei Steshenko <address@hidden>
>> wrote:
>>>  Hello,
>>>
>>>  strictly speaking, it's not an Octave-specific question, but an
>> algorithmic one.
>>>
>>>  Suppose there is a measured function Y(X). In Octave terms X is a vector
>> with N elements.
>>>
>>>  Suppose there are fixed points Xf such that
>>>
>>>  X(1) <= Xf(1)
>>>  Xf(end) <= X(end).
>>>
>>>  The Xf points are more sparse than X.
>>>
>>>
>>>  The Xf points are fixed, i.e. one can't change them as he/she pleases.
>>>
>>>  The goal is to find piecewise-linear function Yf(Xf) which best fits Y(X).
>>>
>>>  I.e. for each two Xf(k), Xf(k+1) pair of points to find a piece of straight
>> line defined by Yf(k), Yf(k+1)pair of points such that the whole Yf fits Y
>> pretty well.
>>>
>>>  Best fitting I'm interested in is according to minimum of sum(abs(Y -
>> Yf_interpolated)). The Yf_interpolated is linear interpolated Yf on X, so
>> dimensions of Y and Yf_interpolated match.
>>>
>>>
>>>  I did some quick web searching and my impression is that there is no
>> universally adopted algorithm for this task, but there is a number solutions,
>> including some for R-language.
>>>
>>>  I myself wrote a straightforward brute force implementation which works
>> pretty well and acceptably fast for me.
>>>
>>>  Anyway, I'm writing this Email in the hope to be educated by the
>> community - maybe there are already more elegant wheels than the one I've
>> invented.
>>>
>>>  Thanks,
>>>    Sergei.
>>>
>>>  _______________________________________________
>>>  Help-octave mailing list
>>> address@hidden
>>>  https://mailman.cae.wisc.edu/listinfo/help-octave
>>
>> Hi,
>>
>> What you describe is also known as Langrange (or linear)
>> interpolation. You can use interp1 with the option linear
>> an example
>>
>> t=linspace(0,2*pi,100);
>> ts=linspace(0,2*pi,10);
>> ys=sin(ts);
>> y=interp1(ts,ys,t,'linear');
>> plot(t,y,'.',ts,ys,'o',t,sin(t),'-')
>>
>> I hope this is what you were asking (to restrict the inteprolation to
>> a subinterval, you could use lookup function before the
>> interpolation).
>>
>>
>>
>> --
>> M. Sc. Juan Pablo Carbajal
>
>
> No, it is _not_ what i am asking, even though I am using a 'interp1' a lot in 
> my code.
>
> This is _not_ interpolation. I am not interested in my "curve" going through 
> certain points, I am interested in _fitting_, i.e. most likely the resulting 
> Yf will _not_ go through any of Y points.
>
> I suggest to enter the subject of this Email into any web search engine - 
> there will be quite a few matches.
>
> Regards,
>   Sergei.

In that case there are at least two options.
- You could formulate your problem as an optimization. Using interp1
to compute the error of you nodes (Xf,Yf) respect to the well-sampled
(X,Y).
- You could solve a least-square-error problem on your data (can be an
optim also, but can by replaced by operations on matrices, which is
faster). Basically you are projecting the curve (X,Y) into the space
defined by the triangular functions (aka tent function
http://en.wikipedia.org/wiki/Tent_function) with centers in Xf. Put
this as a matrix multiplication and use the \ operator or pinv.

What about know? Do I feel you?


-- 
M. Sc. Juan Pablo Carbajal
-----
PhD Student
University of Zürich
http://ailab.ifi.uzh.ch/carbajal/