Paul.Rustomji at csiro.au
2008-Jun-02 04:04 UTC
[R] Missing "spline_coef" DLL and Rob Hyndmans monotonic interpolator
Hello R help I have been trying to use Rob Hyndman's monotonically increasing spline function. But like another user or two seem have a problem with a missing DLL (namely "spline_coef"). None of the previous help postings seemed to have any solutions to this problem. As per a Ripley suggestion I have deleted all previous versions of R and reinstalled R 2.7.0 and the problem persists. Thanks Paul. x <- seq(0,4,l=20) y <- sort(rnorm(20)) plot(x,y) lines(spline(x, y, n = 201), col = 2) # Not necessarily monotonic lines(cm.spline(x, y, n = 201), col = 3)>Error in .C("spline_coef", method = as.integer(method), n = nx, x = x,: C symbol name "spline_coef" not in DLL for package "stats" Cm.spline code from http://www-personal.buseco.monash.edu.au/~hyndman/Rlibrary/interpcode.R ####################################### cm.spline <- function (x, y = NULL, n = 3 * length(x), method = "fmm", xmin = min(x), xmax = max(x), gulim=0) # modification of stats package spline to produce co-monotonic # interpolant by Hyman Filtering. if gulim!=0 then it is taken as the upper # limit on the gradient, and interpolant is gradient limited rather than # monotonic. Modifications from R stats package spline() # (c) Simon N. Wood & Rob J. Hyndman. 2002. { x <- xy.coords(x, y) y <- x$y x <- x$x nx <- length(x) method <- match(method, c("periodic", "natural", "fmm")) if (is.na(method)) stop("spline: invalid interpolation method") dx <- diff(x) if (any(dx < 0)) { o <- order(x) x <- x[o] y <- y[o] } if(any(diff(y)<0)) stop("Data are not monotonic") if (method == 1 && y[1] != y[nx]) { warning("spline: first and last y values differ - using y[1] for both") y[nx] <- y[1] } z <- .C("spline_coef", method = as.integer(method), n = nx, x = x, y = y, b = double(nx), c = double(nx), d = double(nx), e = double(if (method == 1) nx else 0), PACKAGE = "stats") z$y <- z$y-z$x*gulim # trick to impose upper z$b <- z$b-gulim # limit on interpolator gradient z <- hyman.filter(z) # filter gradients for co-monotonicity z$y <- z$y+z$x*gulim # undo trick z$b <- z$b+gulim # transformation z <- spl.coef.conv(z) # force other coefficients to consistency u <- seq(xmin, xmax, length.out = n) .C("spline_eval", z$method, nu = length(u), x = u, y = double(n), z$n, z$x, z$y, z$b, z$c, z$d, PACKAGE = "stats")[c("x","y")] } cm.splinefun <- function(x, y = NULL, method = "fmm",gulim=0) # modification of stats package splinefun to produce co-monotonic #interpolant by Hyman Filtering. if gulim!=0 then it is taken as the upper #limit on the gradient, and intrpolant is gradient limited rather than # monotonic. Modifications from R stats package splinefun() # (c) Simon N. Wood 2002 { x <- xy.coords(x, y) y <- x$y x <- x$x n <- length(x) method <- match(method, c("periodic", "natural", "fmm")) if (is.na(method)) stop("splinefun: invalid interpolation method") if (any(diff(x) < 0)) { z <- order(x) x <- x[z] y <- y[z] } if(any(diff(y)<0)) stop("Data are not monotonic") if (method == 1 && y[1] != y[n]) { warning("first and last y values differ in spline - using y[1] for both") y[n] <- y[1] } z <- .C("spline_coef", method = as.integer(method), n = n, x = as.double(x), y = as.double(y), b = double(n), c double(n), d = double(n), e = double(if (method == 1) n else 0), PACKAGE = "stats") z$y <- z$y-z$x*gulim # trick to impose upper z$b <- z$b-gulim # limit on interpolator gradient z <- hyman.filter(z) # filter gradients for co-monotonicity z$y <- z$y+z$x*gulim # undo trick z$b <- z$b+gulim # transformation z <- spl.coef.conv(z) # force other coefficients to consistency rm(x, y, n, method) function(x) { .C("spline_eval", z$method, length(x), x = as.double(x), y double(length(x)), z$n, z$x, z$y, z$b, z$c, z$d, PACKAGE = "stats")$y } } spl.coef.conv <- function(z) # takes an object z containing equally lengthed arrays # z$x, z$y, z$b, z$c, z$d defining a cubic spline interpolating # z$x, z$y and forces z$c and z$d to be consistent with z$y and # z$b (gradient of spline). This is intended for use in conjunction # with Hyman's monotonicity filter. # Note that R's spline routine has s''(x)/2 as c and s'''(x)/6 as d. # (c) Simon N. Wood { n <- length(z$x) h <- z$x[2:n]-z$x[1:(n-1)] y0 <- z$y[1:(n-1)];y1 <- z$y[2:n] b0 <- z$b[1:(n-1)];b1 <- z$b[2:n] cc <- -(3*(y0-y1)+(2*b0+b1)*h)/h^2 c1 <- (3*(y0[n-1]-y1[n-1])+(b0[n-1]+2*b1[n-1])*h[n-1])/h[n-1]^2 dd <- (2*(y0-y1)/h+b0+b1)/h^2 d1 <- dd[n-1] z$c <- c(cc,c1);z$d <- c(dd,d1) z } hyman.filter <- function(z) # Filters cubic spline function to yield co-monotonicity in accordance # with Hyman (1983) SIAM J. Sci. Stat. Comput. 4(4):645-654, z$x is knot # position z$y is value at knot z$b is gradient at knot. See also # Dougherty, Edelman and Hyman 1989 Mathematics of Computation # 52: 471-494. (c) Simon N. Wood { n <- length(z$x) ss <- (z$y[2:n]-z$y[1:(n-1)])/(z$x[2:n]-z$x[1:(n-1)]) S0 <- c(ss[1],ss) S1 <- c(ss,ss[n-1]) sig <- z$b ind <- (S0*S1>0) sig[ind] <- S1[ind] ind <- (sig>=0) if(sum(ind)) z$b[ind] <- pmin(pmax(0,z$b[ind]),3*pmin(abs(S0[ind]),abs(S1[ind]))) ind <- !ind if(sum(ind)) z$b[ind] <- pmax(pmin(0,z$b[ind]),-3*pmin(abs(S0[ind]),abs(S1[ind]))) z } ####################################### Paul Rustomji Rivers and Estuaries CSIRO Land and Water GPO Box 1666 Canberra ACT 2601 ph +61 2 6246 5810 mobile 0406 375 739