[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: proposal for new m-file function
From: |
Ben Abbott |
Subject: |
Re: proposal for new m-file function |
Date: |
Sat, 24 Mar 2012 10:33:54 -0400 |
On Mar 23, 2012, at 9:23 PM, Martin Helm wrote:
> Am 24.03.2012 01:25, schrieb Ben Abbott:
>
>> On Mar 22, 2012, at 11:16 AM, Ben Abbott wrote:
>>
>>> On Mar 22, 2012, at 10:37 AM, Ben Abbott wrote:
>>>
>>>> On Mar 22, 2012, at 10:12 AM, Martin Helm wrote:
>>>>
>>>>> Am 22.03.2012 13:17, schrieb Ben Abbott:
>>>>>> On Mar 22, 2012, at 6:31 AM, Martin Helm wrote:
>>>>>>
>>>>>>> Am 22.03.2012 10:43, schrieb Carlo de Falco:
>>>>>>>> 2012/3/22<address@hidden
>>>>>>>> <mailto:address@hidden>>
>>>>>>>>
>>>>>>>> Message: 2
>>>>>>>> Date: Wed, 21 Mar 2012 20:33:01 -0400
>>>>>>>> From: Ben Abbott<address@hidden<mailto:address@hidden>>
>>>>>>>> To: octave maintainers mailing list<address@hidden
>>>>>>>> <mailto:address@hidden>>
>>>>>>>> Subject: proposal for new m-file function
>>>>>>>> Message-ID:<address@hidden
>>>>>>>> <mailto:address@hidden>>
>>>>>>>> Content-Type: text/plain; charset="us-ascii"
>>>>>>>>
>>>>>>>> A user's question about fixed points piecewise-linear fitting led
>>>>>>>> to the development of a function looks to be a good fit for
>>>>>>>> Octave's core.
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> https://mailman.cae.wisc.edu/pipermail/help-octave/2012-March/050900.html
>>>>>>>>
>>>>>>>> The idea was to fit a piece-wise polynomial to a set of data.
>>>>>>>> After some discussion between myself and Martin Helm, the attached
>>>>>>>> ppfit.m was produced. The name was chosen to match the existing
>>>>>>>> ppval().
>>>>>>>>
>>>>>>>> The function provides a least-squares fit of a 1D interpolation
>>>>>>>> with specified break positions to a set of data.
>>>>>>>>
>>>>>>>> Demos and tests are included.
>>>>>>>>
>>>>>>>> Any concern about adding this to Octave's core ?
>>>>>>>>
>>>>>>>> Ben
>>>>>>>>
>>>>>>>> googling for the name 'ppfit' I found this function:
>>>>>>>>
>>>>>>>> https://www.assembla.com/code/zaxxon_scripts/subversion/nodes/trunk/project/scripts/ppfit.m
>>>>>>>>
>>>>>>>> which seems to do the same as yours with the additional option of
>>>>>>>> returning a spline with N continuous derivatives
>>>>>>>> Would it be possible to add that option to your code?
>>>>>>>> The license of the linked function looks like BSD so it should be no
>>>>>>>> problem to get the code from there.
>>>>>>>> c.
>>>>>>> Wouldn't it be much cleaner to add that additional fitting options
>>>>>>> (quadratics and higher order splines) to interp1 (without breaking
>>>>>>> matlab compatibility of course)?
>>>>>>> I think that is the place where such functionality should naturally
>>>>>>> live. I could look at it over the weekend and propose a patch for it.
>>>>>> Ok. I'll wait on your changeset.
>>>>>>
>>>>>> Is it possible for interp1 to return the available methods and whether
>>>>>> the underlying interpolants are linear functions ?
>>>>>>
>>>>>> Ben
>>>>> What do you have in mind, is something like interp1("methods") which
>>>>> returns some cellarray for example with the method names and a flag
>>>>> indicating linear/nonlinear an option?
>>>> Something like below?
>>>>
>>>> [methods, linear] = interp1 ("methods")
>>> I just realized I can easily check if the interpolant is linear by ...
>>>
>>> yi = interp1 (xi, eye (numel (xi)), x, "method");
>>> s = sum (yi, 2);
>>> s = s / mean (s);
>>> if (sqrt (mean (abs (s).^2))< sqrt (eps))
>>> islinear = true;
>>> endif
>>>
>>> So, I'll only need to have access to the different methods.
>>>
>>> methods = interp1 ("methods");
>>>
>>> Ben
>> Martin,
>>
>> I've been giving this more thought ... I'm not sure what you had in mind,
>> but after test driving the ppfit.m that Carlo linked to, I'm inclined to;
>>
>> (1) Rename Jonas Lundgren's ppfit.m to __ppfit__.m and then call it from
>> interp1.m and ppfit.m to obtain the piecewise-polynomials for a specified
>> order.
>>
>> (2) Modify interp1.m to allow the method to be specified by "linear",
>> "nearest", "spline", "cubic", "pchip", or numeric value which would be
>> interpreted as the order of polynomial used as the interpolant. Thus, y1 and
>> y2 below would be equivalent.
>>
>> y1 = interp1(xb, yb, x, n);
>>
>> pp = __pfit__ (xb, yb, xb, n);
>> y2 = ppval (pp, x);
>>
>> (3) Modify my ppfit.m to accept the same methods as interp1.m (trivial). I'd
>> also like to support __ppfit__() as it is. Thus, pp1 and pp2 below would be
>> equivalent.
>>
>> pp1 = __ppfit __(x, y, xb, n);
>>
>> pp2 = ppfit (x, y, xb, method, weights, "global");
>>
>> Thoughts ?
>>
>> Ben
>
> I found today some time to test the ppfit function from JonasLundgren myself
> and came more and more to the conclusion that I cannot beat it - not even the
> way it can handle the pure interpolation which is included in his function as
> a subset of its features, while working myself through some algorithms for
> higher order splines.
> The way his implementation handles and constructs the basis for the splines
> is very clever.
> I fully agree that what you propose is the best solution.
> What I had myself in mind does not even come near to it.
I agree. Jonas Lundgren's algorithm is impressive.
I also spent some time with Lundgren's script. I had hoped to find a way to do
the equivalent of the command below in an efficient way.
a = interp1 (xb, eye (numel (xb)), x, method, "extrap");
Rather than modifying his script, perhaps it is best just to do something like
...
yb = eye (numel (xb));
yb = mat2cell (yb, size (yb, 1), ones (1, size (yb, 2)));
func = @(yb) ppval (ppfit (xb, yb, xb, order), x)(:);
a = horzcat (cellfun (func, yb, "UniformOutput", false){:});
In any event, I'll email Jonas and let him know we plan to use his script.
Ben
- Re: proposal for new m-file function, (continued)
- Re: proposal for new m-file function, Carlo de Falco, 2012/03/22
- Re: proposal for new m-file function, Martin Helm, 2012/03/22
- Re: proposal for new m-file function, Ben Abbott, 2012/03/22
- Re: proposal for new m-file function, Martin Helm, 2012/03/22
- Re: proposal for new m-file function, Ben Abbott, 2012/03/22
- Re: proposal for new m-file function, Martin Helm, 2012/03/22
- Re: proposal for new m-file function, Ben Abbott, 2012/03/22
- Re: proposal for new m-file function, Ben Abbott, 2012/03/23
- Re: proposal for new m-file function, Martin Helm, 2012/03/23
- Re: proposal for new m-file function,
Ben Abbott <=
- Re: proposal for new m-file function, Ben Abbott, 2012/03/25
- Re: proposal for new m-file function, Jordi GutiƩrrez Hermoso, 2012/03/26
- Re: proposal for new m-file function, c., 2012/03/26
- Re: proposal for new m-file function, Ben Abbott, 2012/03/26
- Re: proposal for new m-file function, c., 2012/03/26
- Re: proposal for new m-file function, Martin Helm, 2012/03/26
- Re: proposal for new m-file function, Ben Abbott, 2012/03/26
- Re: proposal for new m-file function, Ben Abbott, 2012/03/27
- Re: proposal for new m-file function, Ben Abbott, 2012/03/29