HI there, I have been trying to use a code posted on R help to be able to 
calculate area under the curve for complicated data points and there seems to 
be an issue with the code: no "b" object found. I am not a good R user and 
can't find were the problem is. Any help? Thanks!!
This is the code ( as a test run I gave it this info because I know the answer:
x<-seq(1:50)
 y<-seq(1:50)



simp <- function (y, a = NULL, b = NULL, x = NULL, n = 200)
{
     if (is.null(a) | is.null(b)) {
         if (is.null(x))
             stop("No x values provided to integrate over.\n")
     }
     else {
         x <- c(a, b)
     }
     fff <- 1
     if (length(x) == 2) {
         if (x[1] == x[2])
             return(0)
         if (x[2] < x[1]) {
             fff <- -1
             x <- rev(x)
         }
         x <- seq(x[1], x[2], length = n)
         if (is.function(y))
             y <- y(x)
         else {
             cat("y must be a function when x is\n")
             cat("of length equal to 2.\n")
             stop("Bailing out.\n")
         }
         equisp <- TRUE
     }
     else {
         if (is.function(y))
             y <- y(x)
         else if (length(y) != length(x))
             stop("Mismatch in lengths of x and y.\n")
         s <- order(x)
         x <- x[s]
         ddd <- diff(x)
         if (any(ddd == 0))
             stop("Gridpoints must be distinct.\n")
         equisp <- isTRUE(all.equal(diff(ddd), rep(0, length(ddd) - 1)))        
 y <- y[s]
     }
     n <- length(x) - 1
     if (equisp) {
         old.op <- options(warn = -1)
         on.exit(options(old.op))
         M <- matrix(y, nrow = n + 2, ncol = 4)[1:(n - 2), ]
         h <- x[2] - x[1]
         fc <- h * c(-1, 13, 13, -1)/24
         aa <- apply(t(M) * fc, 2, sum)
         a1 <- h * sum(y[1:3] * c(5, 8, -1))/12
         an <- h * sum(y[(n - 1):(n + 1)] * c(-1, 8, 5))/12
         return(fff * sum(c(a1, aa, an)))
     }
     m <- n%/%2
     i <- 1:(m + 1)
     a <- x[2 * i] - x[2 * i - 1]
     i <- 1:m
     b <- x[2 * i + 1] - x[2 * i]
     o <- (a[i] * b + 2 * a[i] * a[i] - b * b)/(6 * a[i])
     p <- (a[i] + b)^3/(6 * a[i] * b)
     q <- (a[i] * b + 2 * b * b - a[i] * a[i])/(6 * b)
     k <- numeric(n + 1)
     k[1] <- o[1]
     i <- 1:(m - 1)
     k[2 * i] <- p[i]
     k[2 * i + 1] <- q[i] + o[-1]
     if (n > 2 * m) {
         aa <- a[m + 1]
         bb <- b[m]
         den <- 6 * bb * (bb + aa)
         k[2 * m] <- p[m] - (aa^3)/den
         k[2 * m + 1] <- q[m] + (aa^3 + 4 * bb * aa^2 + 3 * aa *
             bb^2)/den
         k[2 * m + 2] <- (2 * bb * aa^2 + 3 * aa * bb^2)/den
     }
     else {
         k[2 * m] <- p[m]
         k[2 * m + 1] <- q[m]
     }
     fff * sum(k * y)
}


--
   Mariana Varela, DVM, PhD
MRC-University of Glasgow Centre for Virus Research
Institute of Infection, Immunity and Inflammation
College of Medical, Veterinary and Life Sciences
Garscube Estate, 464 Bearsden Road
Henry Wellcome Building, room 436
Glasgow, G61 1QH
Scotland (UK)
Phone: +44 (0) 141 330 2196


        [[alternative HTML version deleted]]

______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Reply via email to