Dear all, Can anyone take a look at my program below? There are two functions: f1 (lambda,z,p1) and f2(p1,cl, cu). I fixed p1=0.15 for both functions. For any fixed value of lambda (between 0.01 and 0.99), I solve f1(p1=0.15, lambda=lambda, z)=0 for the corresponding cl and cu values. Then I plug the calculated cl and cu back into the function f2. Eventually, I want to find the lambda value and the corresponding cl and cu values that would make f2=0.1. The result of this program does not seem to match the answer I have. Can some one give me some hint? Thank you very much. Hannah
u1 <- -3 u2 <- 4 f1 <- function(lambda,z,p1){ lambda*(p1*exp(u1*z-u1^2/2)+(0.2-p1)*exp(u2*z-u2^2/2))-(1-lambda)*0.8} f2 <- function(p1,cl, cu){ 0.8*(pnorm(cl)+(1-pnorm(cu)))/(0.8*(pnorm(cl)+(1-pnorm(cu)))+p1*(pnorm(cl- u1)+(1-pnorm(cu-u1)))+(0.2-p1)*(pnorm(cl-u2)+(1-pnorm(cu-u2))))} p1 <- 0.15 lam <- seq(0.01,0.99, by=0.001) x1 <- numeric(length(lam)) for (i in 1:length(lam)){ cl <- uniroot(f1, lower =-10, upper = 0, tol = 1e-10,p1=p1,lambda=lam[i])$root cu <- uniroot(f1, lower =0, upper = 10, tol = 1e-10,p1=p1,lambda=lam[i])$root x1[i]<- f2(p1=p1, cl=cl, cu=cu) } k <- 1 while(k<length(lam) && x1[k]<=0.1){ k=k+1 } k<-k-1;k lower <- uniroot(f1, lower =-10, upper = 0, tol = 1e-10,p1=p1,lambda=lam[k])$root upper <- uniroot(f1, lower =0, upper = 10, tol = 1e-10,p1=p1,lambda=lam[k])$root res <- c(lower, upper) [[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.