[Rd] uniroot violates bounds?
Serguei Sokol
@oko| @end|ng |rom |n@@-tou|ou@e@|r
Mon Feb 20 10:59:22 CET 2023
Le 18/02/2023 à 21:44, J C Nash a écrit :
> I wrote first cut at unirootR for Martin M and he revised and put in
> Rmpfr.
>
> The following extends Ben's example, but adds the unirootR with trace
> output.
>
> c1 <- 4469.822
> c2 <- 572.3413
> f <- function(x) { c1/x - c2/(1-x) }; uniroot(f, c(1e-6, 1))
> uniroot(f, c(1e-6, 1))
> library(Rmpfr)
> unirootR(f, c(1e-6, 1), extendInt="no", trace=1)
>
> This gives more detail on the iterations, and it looks like the Inf is
> the
> issue. But certainly we could do more to avoid "gotchas" like this. If
> someone is interested in some back and forth, I'd be happy to give it a
> try, but I think progress would be better with more than one contributor.
For me, the following fix makes the job :
--- nlm.R.old 2018-09-25 10:44:49.000000000 +0200
+++ nlm.R 2023-02-20 10:46:39.893542531 +0100
@@ -143,14 +143,14 @@
if(check.conv) {
val <- tryCatch(.External2(C_zeroin2, function(arg) f(arg, ...),
- lower, upper, f.lower, f.upper,
+ lower, upper, f.low., f.upp.,
tol, as.integer(maxiter)),
warning = function(w)w)
if(inherits(val, "warning"))
stop("convergence problem in zero finding: ",
conditionMessage(val))
} else {
val <- .External2(C_zeroin2, function(arg) f(arg, ...),
- lower, upper, f.lower, f.upper,
+ lower, upper, f.low., f.upp.,
tol, as.integer(maxiter))
}
iter <- as.integer(val[2L])
Best,
Serguei.
>
> Best,
>
> John Nash
>
> On 2023-02-18 12:28, Ben Bolker wrote:
>> c1 <- 4469.822
>> c2 <- 572.3413
>> f <- function(x) { c1/x - c2/(1-x) }; uniroot(f, c(1e-6, 1))
>> uniroot(f, c(1e-6, 1))
>>
>>
>> provides a root at -6.00e-05, which is outside of the specified
>> bounds. The default value of the "extendInt" argument to uniroot()
>> is "no", as far as I can see ...
>>
>> $root
>> [1] -6.003516e-05
>>
>> $f.root
>> [1] -74453981
>>
>> $iter
>> [1] 1
>>
>> $init.it
>> [1] NA
>>
>> $estim.prec
>> [1] 6.103516e-05
>>
>>
>> I suspect this fails because f(1) (value at the upper bound) is
>> infinite, although setting interval to c(0.01, 1) does work/give a
>> sensible answer ... (works for a lower bound of 1e-4, fails for 1e-5
>> ...)
>>
>> Setting the upper bound < 1 appears to avoid the problem.
>>
>> For what it's worth, the result has an "init.it" component, but the
>> only thing the documentation says about it is " component ‘init.it’
>> was added in R 3.1.0".
>>
>> And, I think (?) that the 'trace' argument only produces any
>> output if the 'extendInt' option is enabled?
>>
>> Inspired by
>> https://stackoverflow.com/questions/75494696/solving-a-system-of-non-linear-equations-with-only-one-unknown/75494955#75494955
>>
>> cheers
>> Ben Bolker
>>
>> ______________________________________________
>> R-devel using r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-devel
>
> ______________________________________________
> R-devel using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
More information about the R-devel
mailing list