From: Terry W. Schulz <terrywschulz_at_comcast.net>

Date: Sun 27 Aug 2006 - 08:42:17 EST

n <- length(count)

mu <- -8.202900475

sigma <- 1.609437912

for(i in c(1:n))

*{
*

integrand <- function(x)

{(x*V[[i]])^c[[i]]*exp(-x*V[[i]])/factorial(c[[i]])*0.3989422804014327/

Date: Sun 27 Aug 2006 - 08:42:17 EST

Hello,

I am a novice R user and am having difficulty retrieving the values from 21 iterations of the R function integrate. The only way I have found is to do a write.table and then a read.table as shown in the code below. I would rather capture the 21 values inside the braces ( sapply might work, but I can't set it up without getting an error in function) so I could compute the negative log-likelihood inside the braces. My goal is to then use optim to find the best fitting mu and sigma. In the code that follows mu and sigma are given, but normally they would not be known.

Any help on y$value capture would be appreciated. I would be ecstatic for help on implementation of the optim function using finite differences for the following code.

By the way in case anyone is interested function(x) is the Poisson-lognormal pdf.

q <- read.table(file="c:/Bayes/FigA3-1.prn",header=TRUE)

attach(q) V <- c(volume) c <- c(count)

n <- length(count)

mu <- -8.202900475

sigma <- 1.609437912

for(i in c(1:n))

integrand <- function(x)

{(x*V[[i]])^c[[i]]*exp(-x*V[[i]])/factorial(c[[i]])*0.3989422804014327/

(x*sigma)*exp(-.5*((log(x)-mu)/sigma)^2)}y <- integrate(integrand, lower = 0, upper = Inf) write.table(y[[1]],file="c:/Bayes/PLNA3-1.out",append=TRUE,

quote=FALSE,row.names=FALSE,col.names=FALSE)
}

z <-read.table(file="c:/Bayes/PLNA3-1.out")
z <- c(z)

LL <- log(c(z$V1))

#negative Log-Likelihood

print(-sum(LL),digits=6)

file.remove("c:/Bayes/PLNA3-1.out")

-- Terry W. Schulz Applied Statistician 1218 Pomegranate Lane Golden, CO 80401 USA terrywschulz@comcast.net (303) 526-1461 ______________________________________________ R-help@stat.math.ethz.ch 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.Received on Sun Aug 27 08:51:18 2006

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.1.8, at Sun 27 Aug 2006 - 10:22:38 EST.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*