uniroot {stats}  R Documentation 
The function uniroot
searches the interval from lower
to upper
for a root (i.e., zero) of the function f
with
respect to its first argument.
Setting extendInt
to a non"no"
string, means searching
for the correct interval = c(lower,upper)
if sign(f(x))
does not satisfy the requirements at the interval end points; see the
‘Details’ section.
uniroot(f, interval, ...,
lower = min(interval), upper = max(interval),
f.lower = f(lower, ...), f.upper = f(upper, ...),
extendInt = c("no", "yes", "downX", "upX"), check.conv = FALSE,
tol = .Machine$double.eps^0.25, maxiter = 1000, trace = 0)
f 
the function for which the root is sought. 
interval 
a vector containing the endpoints of the interval to be searched for the root. 
... 
additional named or unnamed arguments to be passed
to 
lower, upper 
the lower and upper end points of the interval to be searched. 
f.lower, f.upper 
the same as 
extendInt 
character string specifying if the interval

check.conv 
logical indicating whether a convergence warning of the
underlying 
tol 
the desired accuracy (convergence tolerance). 
maxiter 
the maximum number of iterations. 
trace 
integer number; if positive, tracing information is produced. Higher values giving more details. 
Note that arguments after ...
must be matched exactly.
Either interval
or both lower
and upper
must be
specified: the upper endpoint must be strictly larger than the lower
endpoint.
The function values at the endpoints must be of opposite signs (or
zero), for extendInt="no"
, the default. Otherwise, if
extendInt="yes"
, the interval is extended on both sides, in
search of a sign change, i.e., until the search interval [l,u]
satisfies f(l) \cdot f(u) \le 0
.
If it is known how f
changes sign at the root
x_0
, that is, if the function is increasing or decreasing there,
extendInt
can (and typically should) be specified as
"upX"
(for “upward crossing”) or "downX"
,
respectively. Equivalently, define S := \pm 1
, to
require S = \mathrm{sign}(f(x_0 + \epsilon))
at the solution. In that case, the search interval [l,u]
possibly is extended to be such that S\cdot f(l)\le 0
and S \cdot f(u) \ge 0
.
uniroot()
uses Fortran subroutine ‘"zeroin"’ (from Netlib)
based on algorithms given in the reference below. They assume a
continuous function (which then is known to have at least one root in
the interval).
Convergence is declared either if f(x) == 0
or the change in
x
for one step of the algorithm is less than tol
(plus an
allowance for representation error in x
).
If the algorithm does not converge in maxiter
steps, a warning
is printed and the current approximation is returned.
f
will be called as f(x, ...)
for a numeric value
of x.
The argument passed to f
has special semantics and used to be
shared between calls. The function should not copy it.
A list with at least four components: root
and f.root
give the location of the root and the value of the function evaluated
at that point. iter
and estim.prec
give the number of
iterations used and an approximate estimated precision for
root
. (If the root occurs at one of the endpoints, the
estimated precision is NA
.)
Further components may be added in future: component init.it
was added in R 3.1.0.
Based on ‘zeroin.c’ in https://www.netlib.org/c/brent.shar.
Brent, R. (1973) Algorithms for Minimization without Derivatives. Englewood Cliffs, NJ: PrenticeHall.
polyroot
for all complex roots of a polynomial;
optimize
, nlm
.
require(utils) # for str
## some platforms hit zero exactly on the first step:
## if so the estimated precision is 2/3.
f < function (x, a) x  a
str(xmin < uniroot(f, c(0, 1), tol = 0.0001, a = 1/3))
## handheld calculator example: fixed point of cos(.):
uniroot(function(x) cos(x)  x, lower = pi, upper = pi, tol = 1e9)$root
str(uniroot(function(x) x*(x^21) + .5, lower = 2, upper = 2,
tol = 0.0001))
str(uniroot(function(x) x*(x^21) + .5, lower = 2, upper = 2,
tol = 1e10))
## Find the smallest value x for which exp(x) > 0 (numerically):
r < uniroot(function(x) 1e80*exp(x)  1e300, c(1000, 0), tol = 1e15)
str(r, digits.d = 15) # around 745, depending on the platform.
exp(r$root) # = 0, but not for r$root * 0.999...
minexp < r$root * (1  10*.Machine$double.eps)
exp(minexp) # typically denormalized
## uniroot() with new interval extension + checking features: 
f1 < function(x) (121  x^2)/(x^2+1)
f2 < function(x) exp(x)*(x  12)
try(uniroot(f1, c(0,10)))
try(uniroot(f2, c(0, 2)))
##> error: f() .. end points not of opposite sign
## where as 'extendInt="yes"' simply first enlarges the search interval:
u1 < uniroot(f1, c(0,10),extendInt="yes", trace=1)
u2 < uniroot(f2, c(0,2), extendInt="yes", trace=2)
stopifnot(all.equal(u1$root, 11, tolerance = 1e5),
all.equal(u2$root, 12, tolerance = 6e6))
## The *danger* of interval extension:
## No way to find a zero of a positive function, but
## numerically, f(M) becomes zero :
u3 < uniroot(exp, c(0,2), extendInt="yes", trace=TRUE)
## Nonsense example (must give an error):
tools::assertCondition( uniroot(function(x) 1, 0:1, extendInt="yes"),
"error", verbose=TRUE)
## Convergence checking :
sinc < function(x) ifelse(x == 0, 1, sin(x)/x)
curve(sinc, 6,18); abline(h=0,v=0, lty=3, col=adjustcolor("gray", 0.8))
uniroot(sinc, c(0,5), extendInt="yes", maxiter=4) #> "just" a warning
## now with check.conv=TRUE, must signal a convergence error :
uniroot(sinc, c(0,5), extendInt="yes", maxiter=4, check.conv=TRUE)
### Weibull cumulative hazard (example origin, Ravi Varadhan):
cumhaz < function(t, a, b) b * (t/b)^a
froot < function(x, u, a, b) cumhaz(x, a, b)  u
n < 1000
u < log(runif(n))
a < 1/2
b < 1
## Find failure times
ru < sapply(u, function(x)
uniroot(froot, u=x, a=a, b=b, interval= c(1.e14, 1e04),
extendInt="yes")$root)
ru2 < sapply(u, function(x)
uniroot(froot, u=x, a=a, b=b, interval= c(0.01, 10),
extendInt="yes")$root)
stopifnot(all.equal(ru, ru2, tolerance = 6e6))
r1 < uniroot(froot, u= 0.99, a=a, b=b, interval= c(0.01, 10),
extendInt="up")
stopifnot(all.equal(0.99, cumhaz(r1$root, a=a, b=b)))
## An error if 'extendInt' assumes "wrong zerocrossing direction":
uniroot(froot, u= 0.99, a=a, b=b, interval= c(0.1, 10), extendInt="down")