[R] stats:: spline's method could not be monoH.FC

[Martin Maechler]

>>>>> Martin Maechler 
>>>>>     on Mon, 4 May 2020 16:25:02 +0200 writes:

>>>>> Abby Spurdle 
>>>>>     on Sun, 3 May 2020 16:15:17 +1200 writes:

    >>> and just today a colleague asked me about spline interpolation
    >>> with general 2nd derivative boundary conditions
    >>> s''(x_1) = s2_1,  s''(x_n) = s2_2

    > actually I was wrong... I *read* it as the above, but what he
    > really wanted was what the wikipedia page  class "clamped"
    > boundary conditions, i.e., for the *first* derivative

    > s'(x_1) = s1_1,  s'(x_n) = s1_2
   

    >> It should possible via cubic Hermite splines.

    > and indeed that is available via cubic Hermite splines available
    > in R via stats package's  splinefunH()

    > {which I wrote when implementing  "monoH.FC"}

    >> A nontrivial design decision in my package was the computation of
    >> slopes at the endpoints.
    >> (Something which I got wrong, twice...)

    > The "wikiversity" has a very nice small math lecture on this
    > https://en.wikiversity.org/wiki/Cubic_Spline_Interpolation

    > which derives *both* cases (first and 2nd derivative boundary conditions),
    > the only draw back to quickly do it with ('base R') is that I
    > need to "translate" that (2nd derivative values $M_i$) parametrization
    > into either the one used into the (a,b,c,y)-parametrizat of the
    > default spline() / splinefun() methods or the Hilbert spline
    > form for  splinefunH(x[],y[],m[]).

    > Martin

Well, I should have looked a bit further, first : at least the case for

    s''(x_1) = s2_1,    s''(x_n) = s2_2

seems trivially "hidden" in the C code in

R's  src/library/stats/src/splines.c     i.e.
https://svn.r-project.org/R/trunk/src/library/stats/src/splines.c

at the end of natural_spline() we should set c[1] and c[n] to non-zero.
...
and actually I think in spline_eval() for extrapolation (to the
left? and right), there's currently also an assumption that c[1]
 and c[n] are zero.

So as a matter of fact I would ask for patches there (both in C
and in R calling C ... maybe too much for 30 minutes ;-)

Best,
Martin


    >> My guess is that I could write a function in about 60 to 90 minutes,
    >> including all the testing and calculus.

    >> However, I need to note two things:
    >> (1) Cubic Hermite splines do not have a continuous second derivative.

    > well, many don't if you allow any slopes at node. However, if
    > you only set 2 boundary conditions *instead* of the natural
    > spline ones  f''(x_1) = f''(x_1) = 0, 
    > you can still remain in C_2 (i.e. continuous 2nd derivative).

    > (And that's also what the above wikiversity lecture provides).

    >> (2) Specifying the second derivatives (at the endpoints) would prevent
    >> the user from specifying the first derivatives (aka the slopes).

    >> Could you please confirm if the function would still be of interest...?

    > Well, as I know spent enough time reading and thinking, I'd
    > really like to add   method = "clamped" to splinefun() and also
    > the other one where fix the 2nd derivatives (to arbitrary values
    > instead of zero).

    > So if you already have the R code leading up to one of the 2
    > "parametrizations" we use in spline() / splinefun(),  I'd be
    > grateful, if you have the time.

    > Martin