R语言做拟合插值

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

我想用R语言做一个拟合差值。就是我知道变量Ｘ的值然后可以算的变量Y的值。但是x不好由Y算得。所以我想做一个x,y的拟合插值。然后可以由Y的值算得X的值。好像类似于数值分析里面的插值方法。有没有哪位大大帮帮忙呀

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏5 回帖

关键词：R语言 数值分析 有没有

插值算法中比较著名的有牛顿插值，拉格朗日多项式插值，这个是第一型拉格朗日多项式插值算法，不知道是不是你想要的

LagrangePolynomial <- function(x,y) {
  len = length(x)
  if(len != length(y))
    stop("length not equal!")

  if(len < 2)
    stop("dim size must more than 1")

  #pretreat data abd alloc memery
  xx <- paste("(","a -",x,")")
  m <- c(rep(0,len))
  #combin express
  for(i in 1:len) {
    td <- 1
    tm <- "1"
    for(j in 1:len) {
      if(i != j) {
        td <- td*(x - x[j])
        tm <- paste(tm,"*",xx[j])
      }
    }
    tm <- paste(tm,"/",td)
    m<-tm #m <- parse(text=tm)
  }
   
  #combin the exrpession
  m <- paste(m,"*",y)
  r <- paste(m,collapse="+")
   
  #combin the function
  fbody <- paste("{ return(",r,")}")
  f <- function(a) {}
   
  #fill the function's body
  body(f) <- parse(text=fbody)
   
  return(f)
}

算法有两个输入,一个x坐标序列,另一个是y坐标序列
算法返回一个带一个参数的函数
调用次函数得到插值结果

调用方法
>a = 4:6
>b = c(10, 5.25, 1)
>f <- LagrangePolynomial(a,b)
>f(18)
>-11

这个是你要的不

小⑥ 发表于 2014-6-4 16:01
插值算法中比较著名的有牛顿插值，拉格朗日多项式插值，这个是第一型拉格朗日多项式插值算法，不知道是不是 ...

这个我试了一下。差别有点大。我要求出来的数据也是（-1,1）的这个算出来的数差的太大了。不过还是谢谢你啦

那个大的话，你看这个能用不——

require(graphics)

op <- par(mfrow = c(2,1), mgp = c(2,.8,0), mar = 0.1+c(3,3,3,1))
n <- 9
x <- 1:n
y <- rnorm(n)
plot(x, y, main = paste("spline[fun](.) through", n, "points"))
lines(spline(x, y))
lines(spline(x, y, n = 201), col = 2)

y <- (x-6)^2
plot(x, y, main = "spline(.) -- 3 methods")
lines(spline(x, y, n = 201), col = 2)
lines(spline(x, y, n = 201, method = "natural"), col = 3)
lines(spline(x, y, n = 201, method = "periodic"), col = 4)
legend(6, 25, c("fmm","natural","periodic"), col = 2:4, lty = 1)

y <- sin((x-0.5)*pi)
f <- splinefun(x, y)
ls(envir = environment(f))
splinecoef <- get("z", envir = environment(f))
curve(f(x), 1, 10, col = "green", lwd = 1.5)
points(splinecoef, col = "purple", cex = 2)
curve(f(x, deriv = 1), 1, 10, col = 2, lwd = 1.5)
curve(f(x, deriv = 2), 1, 10, col = 2, lwd = 1.5, n = 401)
curve(f(x, deriv = 3), 1, 10, col = 2, lwd = 1.5, n = 401)
par(op)

## Manual spline evaluation --- demo the coefficients :
.x <- splinecoef$x
u <- seq(3, 6, by = 0.25)
(ii <- findInterval(u, .x))
dx <- u - .x[ii]
f.u <- with(splinecoef,
            y[ii] + dx*(b[ii] + dx*(c[ii] + dx* d[ii])))
stopifnot(all.equal(f(u), f.u))

## An example with ties (non-unique  x values):
set.seed(1); x <- round(rnorm(30), 1); y <- sin(pi * x) + rnorm(30)/10
plot(x, y, main = "spline(x,y)  when x has ties")
lines(spline(x, y, n = 201), col = 2)
## visualizes the non-unique ones:
tx <- table(x); mx <- as.numeric(names(tx[tx > 1]))
ry <- matrix(unlist(tapply(y, match(x, mx), range, simplify = FALSE)),
             ncol = 2, byrow = TRUE)
segments(mx, ry[, 1], mx, ry[, 2], col = "blue", lwd = 2)

## An example of monotone interpolation
n <- 20
set.seed(11)
x. <- sort(runif(n)) ; y. <- cumsum(abs(rnorm(n)))
plot(x., y.)
curve(splinefun(x., y.)(x), add = TRUE, col = 2, n = 1001)
curve(splinefun(x., y., method = "monoH.FC")(x), add = TRUE, col = 3, n = 1001)
curve(splinefun(x., y., method = "hyman")   (x), add = TRUE, col = 4, n = 1001)
legend("topleft",
       paste0("splinefun( \"", c("fmm", "monoH.FC", "hyman"), "\" )"),
       col = 2:4, lty = 1, bty = "n")

## and one from Fritsch and Carlson (1980), Dougherty et al (1989)
x. <- c(7.09, 8.09, 8.19, 8.7, 9.2, 10, 12, 15, 20)
f <- c(0, 2.76429e-5, 4.37498e-2, 0.169183, 0.469428, 0.943740,
       0.998636, 0.999919, 0.999994)
s0 <- splinefun(x., f)
s1 <- splinefun(x., f, method = "monoH.FC")
s2 <- splinefun(x., f, method = "hyman")
plot(x., f, ylim = c(-0.2, 1.2))
curve(s0(x), add = TRUE, col = 2, n = 1001) -> m0
curve(s1(x), add = TRUE, col = 3, n = 1001)
curve(s2(x), add = TRUE, col = 4, n = 1001)
legend("right",
       paste0("splinefun( \"", c("fmm", "monoH.FC", "hyman"), "\" )"),
       col = 2:4, lty = 1, bty = "n")

## they seem identical, but are not quite:
xx <- m0$x
plot(xx, s1(xx) - s2(xx), type = "l",  col = 2, lwd = 2,
     main = "Difference   monoH.FC - hyman"); abline(h = 0, lty = 3)

x <- xx[xx < 10.2] ## full range: x <- xx .. does not show enough
ccol <- adjustcolor(2:4, 0.8)
matplot(x, cbind(s0(x, deriv = 2), s1(x, deriv = 2), s2(x, deriv = 2))^2,
        lwd = 2, col = ccol, type = "l", ylab = quote({{f*second}(x)}^2),
        main = expression({{f*second}(x)}^2 ~" for the three 'splines'"))
legend("topright",
       paste0("splinefun( \"", c("fmm", "monoH.FC", "hyman"), "\" )"),
       lwd = 2, col  =  ccol, lty = 1:3, bty = "n")
## --> "hyman" has slightly smaller  Integral f''(x)^2 dx  than "FC",   
这里好几个例子，可以选一个不...呃，总有一款适合你吧...

DM小菜鸟 发表于 2014-12-28 22:20
那个大的话，你看这个能用不——

require(graphics)

谢谢大大

小⑥ 发表于 2014-6-4 16:01
插值算法中比较著名的有牛顿插值，拉格朗日多项式插值，这个是第一型拉格朗日多项式插值算法，不知道是不是 ...

请问有克里金插值的程序么？求教！

DM小菜鸟 发表于 2014-12-28 22:20
那个大的话，你看这个能用不——

require(graphics)

请问有克里金插值的程序么？谢谢！

小⑥ 发表于 2014-6-4 16:01
插值算法中比较著名的有牛顿插值，拉格朗日多项式插值，这个是第一型拉格朗日多项式插值算法，不知道是不是 ...

这个代码运行完会返回Inf，不知道是什么原因

我的用第一个得到的也是Inf，不知道为啥

我也是，你懂了吗