[R] optimize simultaneously two binomials inequalities using nlm( ) or optim( )

Emir Toktar emir.toktar at gmail.com
Wed Aug 6 19:05:32 CEST 2008


Thanks Mr. Graves for your support.

>       I saw your post on 7/29, and I have not seen a reply, 
> so I will attempt a response to the question at the start of 
> your email:  
> 
> obtain the smallest value of 'n' (sample size) satisfying 
> both inequalities: 
> 
> 		    (1-alpha) <= pbinom(c, n, p1) && pbinom(c, 
> n, p2) <= beta, 
> 
> where alpha, p1, p2, and beta are given, and I assume that 
> 'c' is also given, though that's not perfectly clear to me.  

Ok, the 'c' starts from 0 and the limit was sptipulate to lot size limit
(N). When the c value is increased, the OC curve is shifted to the right.
I'm using the lot size concept because the hypergeometric accetance model
that I'm using too.
 
> 
> 	  Since 'n' is an integer, standard optimizers like 
> optim and nlm are not really appropriate.  This sounds like 
> integer programming.  RSiteSearch('integer programming', 
> 'fun') just produced 138 hits for me.  You might find 
> something useful there.  

Thanks, I'm searching...
> 
> 	  However, before I tried that, I'd try simpler things 
> first.  Consider for example the following:  
> 
> c. <- 5
> # I used 'c.' not 'c', because 'c' is the name of a function in R.  
> alpha <- .05
> beta <- .8
> p1 <- .2
> p2 <- .9 
> 
> n <- c.:20
> p.1 <- pbinom(c., n, p1)
> p.2 <- pbinom(c., n, p2)
> 
> good <- (((1-alpha) <= p.1) & (p.2 <= beta))
> min(n[good]) 
> 
> op <- par(mfrow=c(2, 1))
> plot(n, p.1)
> abline(h=1-alpha)
> plot(n, p.2)
> abline(h=beta)
> par(op)
> 
> 	  In this case, n = 6 satisfies both inequalities, 
> though n = 15 does not.  If min(n[good]) = Inf either no 
> solution exists or you need to increase the upper bound of 
> your search range for 'n'.  

I did an interactive method before, but not using an R optimizer method.

findN <- function(N, p2, beta, c.)
{ for (ss in c.:N)	## 'j' is the sample size value
     if (pbinom(c., ss, p2) <= beta) return (ss);    
     return (N); # bad search, limit is the lot size N
} 
findC <- function(N, p1, Alpha, c., k)
{ for (d in c.:N)#'i' is the quantile value (defective)
     if (pbinom(d, k, p1) >= Alpha) return(d);
     return(N); # bad search, limit is the lot size N
} 
recur.binom <- function(N, p1, p2, alpha, beta)
{
 	c. <- 0;	k <- 0;  conf.level <- (1 - alpha)
	while (1) {
	    k <- findN(N, p2, beta, c.);
	    if (pbinom(c., k, p1) >= conf.level)   break
	    else    c. <- findC(N, p1, conf.level, c., k);   }

  cat(paste(" Acceptance Number 'c' : ", c. ," Sample number 'k' : ",
k,"\n"));
  cat(paste(" Confidence Level=",pbinom(c., k, p1) * 100,"%\n"));
  prob <- seq(0.0001,0.15,by=0.001);
  Pb <- pbinom(c., k, prob);
  plot(prob, Pb, type='l');
}
recur.binom(N = 600, p1 = 0.01, p2 = 0.08, alpha = 0.05, beta = 0.05)
 Sample number 'k' :  77  Acceptance Number 'c' :  2 
 Confidence Level= 95.76 %

I'll trying modified this case and use your suggestion, perhaps changing the
p2 definition(Consumer option) as the maximum value and not fixed value used
in my example, the same to alpha and beta values, giving a more flexible
model to optimize.

> 
> 	  If you'd like more help, PLEASE do read the posting 
> guide 'http://www.R-project.org/posting-guide.html' and 
> provide commented, minimal, self-contained, reproducible code.
> 
> 	  Hope this helps.  
> 	  Spencer Graves
> 

Ok, this was my first post. My apologise for the inconvenient. 
Thank you very much!!

Etoktar
+55 (21) 9112-1847
emir.toktar at computer.org
emir.toktar at gmail.com



More information about the R-help mailing list