This task view on numerical mathematics lists R packages and functions
that are useful for solving numerical problems in linear algebra and
analysis. It shows that R is a viable computing environment for
implementing and applying numerical methods, also outside the realm of
statistics.
The task view will
not
cover differential equations,
optimization problems and solvers, or packages and functions operating
on times series, because all these topics are treated extensively in
the corresponding task views
DifferentialEquations,
Optimization, and
TimeSeries. All these task
views together will provide a good selection of what is available in R for
the area of numerical mathematics.
The task view has been created to provide an overview of the topic.
If some packages are missing or certain topics in numerical math
should be treated in more detail, please let the maintainer know.
Numerical Linear Algebra
As statistics is based to a large extent on linear algebra, many
numerical linear algebra routines are present in R, and some only
implicitly. Examples of explicitly available functions are vector and
matrix operations, matrix (QR) decompositions, solving linear equations,
eigenvalues/-vectors, singular value decomposition, or least-squares
approximation.
-
The recommended package
Matrix
provides classes and methods
for dense and sparse matrices and operations on them, for example
Cholesky and Schur decomposition, matrix exponential, or norms and
conditional numbers for sparse matrices.
-
Recommended package
MASS
adds generalized (Penrose)
inverses and null spaces of matrices.
-
expm
computes the exponential, logarithm, and square root
of square matrices, but also powers of matrices or the Frechet
derivative.
expm()
is to be preferred to the function
with the same name in
Matrix.
-
SparseM
provides classes and methods for sparse matrices
and for solving linear and least-squares problems in sparse linear
algebra
-
svd
provides R bindings to state-of-the-art implementations
of singular value decomposition (SVD) and eigenvalue/eigenvector
computations. Package
ssvd
will obtain sparse SVDs using an
iterative thresholding method, while
irlba
will compute
approximate singular values/vectors of large matrices.
-
The packages
geigen
and
QZ
compute generalized
eigenvalues and -vectors for pairs of matrices, and QZ (generalized
Schur) decompositions.
-
eigeninv
generates matrices with a given set of
eigenvalues ('inverse eigenvalue problem').
-
Package
onion
contains routines for manipulating
quaternions and octonians (normed division algebras over the real
numbers); quaternions can be useful for handling rotations in
three-dimensional space.
-
Packages
RcppArmadillo
and
RcppEigen
enable the
integration of the C++ template libraries 'Armadillo' resp. 'Eigen'
for linear algebra applications written in C++ and integrated in R
using
Rcpp
for performance and ease of use.
Special Functions
Many special mathematical functions are present in R, especially
logarithms and exponentials, trigonometric and hyperbolic functions, or
Bessel and Gamma functions. Many more special functions are available in
contributed packages.
-
Package
gsl
provides an interface to the 'GNU Scientific
Library' that contains implementations of many special functions, for
example the Airy and Bessel functions, elliptic and exponential
integrals, the hypergeometric function, Lambert's W function, and
many more.
-
Airy and Bessel functions, for real and complex numbers, are also
computed in package
Bessel, with approximations for large
arguments.
-
Package
pracma
includes special functions, such as
error functions and inverses, incomplete and complex gamma function,
exponential and logarithmic integrals, Fresnel integrals, the
polygamma and the Riemann zeta function.
-
appell
computes Gauss' 2F1 and Appell's F1 hypergeometric
functions for complex parameters and arguments quite accurately.
-
The hypergeometric (and generalized hypergeometric) function, is
computed in
hypergeo, including transformation formulas
and special values of the parameters.
-
Elliptic and modular functions are available in package
elliptic, including the Weierstrass P function and Jacobi's
theta functions.
There are tools for visualizing complex functions.
Polynomials
Function polyroot() in base R determines all zeros of a polynomial,
based on the Jenkins-Traub algorithm. Linear regression function lm()
can perform polynomial fitting when using
poly()
in the model
formula (with option
raw = TRUE).
-
Packages
polynom
and
PolynomF
provide similar
functionality for manipulating univariate polynomials, like
evaluating polynomials (Horner scheme), differentiating or integrating
them, or solving polynomials, i.e. finding all roots (based on an
eigenvalue computation).
-
Package
MonoPoly
fits univariate polynomials to given data,
applying different algorithms.
-
For multivariate polynomials, package
multipol
provides
various tools to manipulate and combine these polynomials of several
variables.
-
Package
mpoly
facilitates symbolic manipulations on
multivariate polynomials, including basic differential calculus
operations on polynomials, plus some Groebner basis calculations.
-
Package
orthopolynom
consists of a collection of functions
to construct orthogonal polynomials and their recurrence relations,
among them Chebyshev, Hermite, and Legendre polynomials, as well as
spherical and ultraspherical polynomials. There are functions to
operate on these polynomials.
Differentiation and Integration
D()
and
deriv()
in base R compute derivatives
of simple expressions symbolically.
Function
integrate()
implements an approach for numerically
integrating univariate functions in R. It applies adaptive Gauss-Kronrod
quadrature and can handle singularities and unbounded domains to a certain
extent.
-
numDeriv
sets the standard for numerical differentiation
in R, providing numerical gradients, Jacobians, and Hessians, computed
by simple finite differences, Richardson extrapolation, or the highly
accurate complex step approach.
-
features
extracts features from functional data, such as
first and second derivatives, or curvature at critical points.
-
Package
pracma
contains functions for computing numerical
derivatives, including Richardson extrapolation or complex step.
fderiv()
computes numerical derivatives of higher orders.
pracma
has several routines for numerical integration:
adaptive Lobatto quadrature, Romberg integration, Newton-Cotes
formulas, Clenshaw-Curtis quadrature rules.
integral2()
integrates functions in two dimensions, also for domains characterized
by polar coordinates or with variable interval limits.
-
Package
gaussquad
contains a collection of functions to
perform Gaussian quadrature, among them Chebyshev, Hermite, Laguerre,
and Legendre quadrature rules, explicitly returning nodes and weights
in each case. Function
gaussquad()
in package
statmod
does a similar job.
-
Adaptive multivariate integration over hyper-rectangles in
n-dimensional space is available in package
cubature
as
function
adaptIntegrate(), based on C code by Stephen
Johnson. The integrand functions can even be multi-valued.
-
Package
SphericalCubature
employs
cubature
to
integrate functions over unit spheres and balls in n-dimensional
space; conversion routines for multivariate, spherical, and polar
coordinates are provided.
-
Multi-dimensional numerical integration is also covered in package
R2Cuba, employing a C library by Thomas Hahn.
With
vegas()
it includes an approach to Monte Carlo
integration based on importance sampling.
-
Package
SparseGrid
provides another approach to
multivariate integration in high-dimensional spaces. It creates sparse
n-dimensional grids that can be used as with quadrature rules.
-
Package
polyCub
holds some routines for numerical
integration over polygonal domains in two dimensions.
Interpolation and Approximation
Base R provides functions
approx()
for constant and linear
interpolation, and
spline()
for cubic (Hermite) spline
interpolation, while
smooth.spline()
performs cubic spline
approximation. Base package splines creates periodic interpolation
splines in function
periodicSpline().
-
Interpolation of irregularly spaced data is possible with the
akima
package:
aspline()
for univariate data,
bicubic()
or
interp()
for data on a 2D
rectangular domain. (This package is distributed under ACM license and
not available for commercial use.)
-
Package
signal
contains several
filters
to smooth
discrete data, notably
interp1()
for linear, spline, and
cubic interpolation,
pchip()
for piecewise cubic Hermite
interpolation, and
sgolay()
for Savitzky-Golay
smoothing.
-
Package
pracma
provides barycentric Lagrange interpolation
(in 1 and 2 dimensions) in
barylag()
resp.
barylag2d(), 1-dim. akima in
akimaInterp(),
and interpolation and approximation of data with rational functions,
i.e. in the presence of singularities, in
ratinterp()
and
rationalfit().
-
sinterp()
in package
stinepack
realizes
interpolation based on piecewise rational functions by applying
Stineman's algorithm. The interpolating function will be monoton in
regions where the specified points change monotonically.
Root Finding
uniroot(), implementing the Brent-Decker algorithm, is the
basic routine in R to find roots of univariate functions. There are
implementations of the bisection algorithm in several contributed
packages. For root finding with higher precision there is function
unirootR()
in the multi-precision package
Rmpfr.
And for finding roots of multivariate functions see the following two
packages:
-
For solving nonlinear systems of equations the
BB
package
provides (non-monotone) Barzilai-Borwein spectral methods in
sane(), including a derivative-free variant in
dfsane(), and multi-start features with sensitivity
analysis.
-
Package
nleqslv
solves nonlinear systems of equations
using alternatively the Broyden or Newton method, supported by
strategies such as line searches or trust regions.
-
ktsolve
defines a common interface for solving a set of
equations with
BB
or
nleqslv.
Discrete Mathematics and Number Theory
Not so many functions are available for computational number theory.
Note that integers in double precision can be represented exactly up to
2^53 - 1, above that limit a multi-precision package such as
gmp
is needed, see below.
-
Package
numbers
provides functions for factorization, prime
numbers, twin primes, primitive roots, modular inverses, extended GCD,
etc. Included are some number-theoretic functions like divisor
functions or Euler's Phi function.
-
contfrac
contains various utilities for evaluating
continued fractions and partial convergents.
-
The
partitions
package enumerates additive partitions
of integers, including restricted and unequal partitions.
-
Package
combinat
generates all permutations or all
combinations of a certain length of a set of elements (i.e. a vector);
it also computes binomial coefficients.
pcenum
provides
methods to
iterate
through all permutations or combinations
without retrieving them all at one time.
-
magic
creates and investigates magical squares.
Multi-Precision Arithmetic and Symbolic Mathematics
-
Multiple precision arithmetic is available in R through package
gmp
that interfaces to the GMP C library. Examples are
factorization of integers, a probabilistic prime number test, or
operations on big rationals -- for which linear systems of equations
can be solved.
-
Multiple precision floating point operations and functions are
provided through package
Rmpfr
using the MPFR and GMP
libraries. Special numbers and some special functions are included,
as well as routines for root finding, integration, and optimization
in arbitrary precision.
-
Package
rSymPy
accesses the symbolic algebra system 'SymPy'
(written in Python) from R.
It supports arbitrary precision computations, linear algebra and
calculus, solving equations, discrete mathematics, and much more.
-
Package
Ryacas
interfaces the computer algebra system
'Yacas'. It supports symbolic and arbitrary precision computations
in calculus and linear algebra.
MATLAB, Octave, and Python Interfaces
Interfaces to numerical computation software such as MATLAB
(commercial) or Octave (free) will be important when solving difficult
numerical problems. Please note that the commercial programs SAS and
Mathematica do have facilities to call R functions.
-
The
matlab
emulation package contains about 30 simple
functions, replicating MATLAB functions, using the respective MATLAB
names and being implemented in pure R.
-
Package
R.matlab
provides tools to read and write MAT
files, which is the MATLAB data format. It also enables a
one-directional interface with a MATLAB process, sending and
retrieving objects through a TCP/IP connection.
-
Package
RcppOctave
provides an interface to Octave (a
MATLAB clone). It enables calling Octave functions, passing variables
between R and Octave, and browsing Octave documentation.
-
Python, through its modules 'NumPy', 'SciPy', and 'pandas', has
elaborate and efficient numerical tools available. R package
rPython
permits calls from R to Python, while project
RPy
(with Python
module 'rpy2') interfaces R from Python.