1

Finding root for Hermite spline I have a question regarding one discussion that already in StackOverflow
get x-value given y-value: general root finding for linear / non-linear interpolation function

They told that for cubic interpolation splines returned by stats::splinefun with methods "fmm", "natrual", "periodic" and "hyman", the function RootSpline3 provides a stable numerical solution. Now if I use "monoH.FC", would the RootSpline3 function work? I have actually tried it, but seems like it doesn't work. Can you tell what is wrong in my code(why the length argument is invalid)? Either my code is wrong or it will not work for this particular method? If yes, what I need to do?).

kne<-c(10,15,18,18,15,14,13,13,15,21,26,39,52,64,70,66,57,40,22,11)
t<-seq(0,1,len=20)
s <- splinefun(t, kne, method = "monoH.FC")

RootSpline3 <- function (s, y0 = 0, verbose = TRUE) {
  ## extract piecewise construction info
  info <- environment(s)$z
  print(info)
  n_pieces <- info$n - 1L
  x <- info$x; y <- info$y
  print(x)
  b <- info$b; c <- info$c; d <- info$d
  ## list of roots on each piece
  xr <- vector("list", n_pieces)
  ## loop through pieces
  i <- 1L
  while (i <= n_pieces) {
    ## complex roots
    croots <- polyroot(c(y[i] - y0, b[i], c[i], d[i]))
    ## real roots (be careful when testing 0 for floating point numbers)
    rroots <- Re(croots)[round(Im(croots), 10) == 0]
    ## the parametrization is for (x - x[i]), so need to shift the roots
    rroots <- rroots + x[i]
    ## real roots in (x[i], x[i + 1])
    xr[[i]] <- rroots[(rroots >= x[i]) & (rroots <= x[i + 1])]
    ## next piece
    i <- i + 1L
  }
  ## collapse list to atomic vector
  xr <- unlist(xr)
  ## make a plot?
  if (verbose) {
    curve(f, from = x[1], to = x[n_pieces + 1], xlab = "x", ylab = "f(x)")
    abline(h = y0, lty = 2)
    points(xr, rep.int(y0, length(xr)))
  }
  ## return roots
  xr
}

RootSpline3(s, 10)
Antora
  • 11
  • 4
  • Please provide a [mcve]. For example, `Data` and `kne` are not defined. And is it possible you have some `n_pieces` which are `NA`? – andrew_reece May 03 '19 at 02:37
  • @andrew_reece, I have made a simple random example to edit the question and as far as I think n_pieces should be a number. Though if you try to print it, it gives NA, but I don't understand why. – Antora May 03 '19 at 02:50
  • @Antora Thanks for playing with my function. Yes, I need to think about this. – Zheyuan Li Jul 26 '22 at 02:51

1 Answers1

0

Depending on which spline method you use, you'll get a different environment. In the case of method = 'monoH.FC', there are no environment variables available (like $z, $c, $d, etc) as there are in the default spline method ('fmm').

However, there is still an x and y argument that can be called with a monoH.FC splinefun(), both of which will produce vectors, and which can be used to plot the spline fit.

Using your data:

s(kne)
# -1870 -14201  -3542  -3751  -2915  -2706  -2497  -3333  -2288  -4169  -5214 
# -7931 -10648 -13156 -14410 -13574 -11693  -8140  -4378  -2079
s(t)
# 10 69 18 19 15 14 13 17 12 21 26 39 52 64 70 66 57 40 22 11

plot(t, kne)
curve(s(x), add = TRUE)

spline plot

andrew_reece
  • 20,390
  • 3
  • 33
  • 58
  • 1
    actually my question is how to find x for y(function value) say 10? I have tried "x[which(s(x) == 10)]", but it gives NULL. – Antora May 03 '19 at 03:51
  • 1
    I mean It gives values for the points we already have, but not for other points(it actually totally makes sense, because we are not doing inverse interpolating). I am just wondering any way to do inverse interpolation for this particular method. I am new in R and couldn't write my own function for inverse interpolation. – Antora May 03 '19 at 04:27