[R] using optimize with two unknowns, e.g. to parameterize a distribution with given confidence interval
Ravi Varadhan
rvaradhan at jhmi.edu
Tue Oct 19 17:09:11 CEST 2010
You cannot use `optimize' when there are two or more parameters to be
optimized. I don’t know if other have suggested any solution to this, but
here are 2 approaches:
# Estimating LCL and UCL separately using `optimize'.
prior.lcl <- function(x, alpha, mean, var) {
a <- abs(plnorm(x, mean, var) - (alpha/2))
return(a)
}
prior.ucl <- function(x, alpha, mean, var) {
b <- abs(plnorm(x, mean, var) - (1-alpha/2))
return(b)
}
optimize(f=prior.lcl, mean=0, var=1, alpha=0.05, interval=c(0,20))
optimize(f=prior.ucl, mean=0, var=1, alpha=0.05, interval=c(0,20))
# Combining LCL and UCL estimation
# Using `optimx' package to illustrate how to solve this using different
optimizers
# This also shows that this problem is not so easy to solve
#
prior <- function(x, alpha, mean, var) {
a <- abs(plnorm(x[1], mean, var) - (alpha/2))
b <- abs(plnorm(x[2], mean, var) - (1-alpha/2))
return(a+b)
}
require(optimx)
optimx(par=c(0,10), fn=prior, mean=0, var=1, alpha=0.05, method=c("Nelder",
"BFGS", "CG", "spg", "bobyqa", "Rvmmin", "nlminb", "Rcgmin", "ucminf"))
optimx(par=c(1,10), fn=prior, mean=0, var=1, alpha=0.05, method=c("Nelder",
"BFGS", "CG", "spg", "bobyqa", "Rvmmin", "nlminb", "Rcgmin", "ucminf"))
Hope this helps,
Ravi.
-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
Behalf Of David LeBauer
Sent: Friday, October 15, 2010 5:31 PM
To: r-help
Subject: [R] using optimize with two unknowns, e.g. to parameterize a
distribution with given confidence interval
Hi,
I would like to write a function that finds parameters of a log-normal
distribution with a 1-alpha CI of (x_lcl, x_ucl):
However, I don't know how to optimize for the two unknown parameters.
Here is my unsuccessful attempt to find a lognormal distribution with
a 90%CI of 1,20:
prior <- function(x_lcl, x_ucl, alpha, mean, var) {
a <- (plnorm(x_lcl, mean, var) - (alpha/2))^2
b <- (plnorm(x_ucl, mean, var) - (1-alpha/2))^2
return(a+b)
}
optimize(fn=prior, interval = c(-5, 100), 1, 20)
I understand that this problem has a closed form solution, but I would
like to make this a general function.
Thanks,
David
--
David LeBauer, PhD
Energy Biosciences Institute
University of Illinois Urbana-Champaign
1206 W. Gregory Drive
Urbana, IL 61801, U.S.A.
______________________________________________
R-help at 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.
More information about the R-help
mailing list