https://github.com/cran/pracma
Raw File
Tip revision: 708a2ad382a163d1eef5af0665e3ae2aad200ced authored by HwB on 21 March 2013, 00:00:00 UTC
version 1.4.5
Tip revision: 708a2ad
lsqnonlin.Rd
\name{lsqnonlin}
\alias{lsqnonlin}
\alias{lsqnonneg}
\alias{lsqcurvefit}
\alias{lsqsep}
\title{
  Nonlinear Least-Squares Fitting
}
\description{
  \code{lsqnonlin} solves nonlinear least-squares problems, including
  nonlinear data-fitting problems, through the Levenberg-Marquardt approach.

  \code{lsqnonneg} solve nonnegative least-squares constraints problem.
}
\usage{
lsqnonlin(fun, x0, options = list(), ...)
lsqnonneg(C, d)

lsqsep(flist, p0, xdata, ydata, const = TRUE)
lsqcurvefit(fun, p0, xdata, ydata)
}
\arguments{
  \item{fun}{User-defined, vector-valued function.}
  \item{x0}{starting point.}
  \item{...}{additional parameters passed to the function.}
  \item{options}{list of options, for details see below.}
  \item{C, d}{matrix and vector such that \code{C x - d} will be
              minimized with \code{x >= 0}.}

  \item{flist}{list of (nonlinear) functions, depending on one extra parameter.}
  \item{p0}{starting parameters.}
  \item{xdata, ydata}{data points to be fitted.}
  \item{const}{logical; shall a constant term be included.}
}
\details{
  \code{lsqnonlin} computes the sum-of-squares of the vector-valued function
  \code{fun}, that is if \eqn{f(x) = (f_1(x), \ldots ,f_n(x))} then
  \deqn{min || f(x) ||_2^2 = min(f_1(x)^2 + \ldots + f_n(x)^2)}
  will be minimized.

  \code{x=lsqnonlin(fun,x0)} starts at point \code{x0} and finds a minimum
  of the sum of squares of the functions described in fun. \code{fun} shall
  return a vector of values and not the sum of squares of the values.
  (The algorithm implicitly sums and squares fun(x).)

  \code{options} is a list with the following components and defaults:
  \itemize{
    \item \code{tau}: used as starting value for Marquardt parameter.
    \item \code{tolx}: stopping parameter for step length.
    \item \code{tolg}: stopping parameter for gradient.
    \item \code{maxeval} the maximum number of function evaluations.
  }
  Typical values for \code{tau} are from \code{1e-6...1e-3...1} with small
  values for good starting points and larger values for not so good or known
  bad starting points.

  \code{lsqnonneg} solves the linear least-squares problem \code{C x - d},
  \code{x} nonnegative, transforming it with the `trick' \code{x --> exp(x)}
  into a nonlinear one and solving it applying \code{lsqnonlin}.

  \code{lsqsep} solves the separable least-squares fitting problem

  \code{y = a0 + a1*f1(b1, x) + ... + an*fn(bn, x)}

  where \code{fi} are nonlinear functions each depending on a single extra
  paramater \code{bi}, and \code{ai} are additional linear parameters that
  can be separated out to solve a nonlinear problem in the \code{bi} alone.

  \code{lsqcurvefit} is simply an application of \code{lsqnonlin} to fitting
  data points. \code{fun(p, x)} must be a function of two groups of variables
  such that \code{p} will be varied to minimize the least squares sum, see
  the example below.
}
\value{
  \code{lsqnonlin} returns a list with the following elements:
  \itemize{
    \item \code{x}: the point with least sum of squares value.
    \item \code{ssq}: the sum of squares.
    \item \code{ng}: norm of last gradient.
    \item \code{nh}: norm of last step used.
    \item \code{mu}: damping parameter of Levenberg-Marquardt.
    \item \code{neval}: number of function evaluations.
    \item \code{errno}: error number, corresponds to error message.
    \item \code{errmess}: error message, i.e. reason for stopping.
  }

  \code{lsqsep} will return the coefficients sparately, \code{a0} for the
  constant term (being 0 if \code{const=FALSE}) and the vectors \code{a} and
  \code{b} for the linear and nonlinear terms, respectively.
}
\note{
  The refined approach, Fletcher's version of the Levenberg-Marquardt
  algorithm, may be added at a later time; see the references.
}
\references{
  Madsen, K., and H. B.Nielsen (2010). Introduction to Optimization and
  Data Fitting. Technical University of Denmark, Intitute of Computer
  Science and Mathematical Modelling.

  Fletcher, R., (1971). A Modified Marquardt Subroutine for Nonlinear Least
  Squares. Report AERE-R 6799, Harwell.
}
\seealso{
  \code{\link{nlm}}, \code{\link{nls}}
}
\examples{
##  Rosenberg function as least-squares problem
x0  <- c(0, 0)
fun <- function(x) c(10*(x[2]-x[1]^2), 1-x[1])
lsqnonlin(fun, x0)

##  Example from R-help
y <- c(5.5199668,  1.5234525,  3.3557000,  6.7211704,  7.4237955,  1.9703127,
       4.3939336, -1.4380091,  3.2650180,  3.5760906,  0.2947972,  1.0569417)
x <- c(1,   0,   0,   4,   3,   5,  12,  10,  12, 100, 100, 100)
# Define target function as difference
f <- function(b)
     b[1] * (exp((b[2] - x)/b[3]) * (1/b[3]))/(1 + exp((b[2] - x)/b[3]))^2 - y
x0 <- c(21.16322, 8.83669, 2.957765)
lsqnonlin(f, x0)        # ssq 50.50144 at c(36.133144, 2.572373, 1.079811)

# nls() will break down
# nls(Y ~ a*(exp((b-X)/c)*(1/c))/(1 + exp((b-X)/c))^2,
#     start=list(a=21.16322, b=8.83669, c=2.957765), algorithm = "plinear")
# Error: step factor 0.000488281 reduced below 'minFactor' of 0.000976563

##  Example: Hougon function
x1 <- c(470, 285, 470, 470, 470, 100, 100, 470, 100, 100, 100, 285, 285)
x2 <- c(300,  80, 300,  80,  80, 190,  80, 190, 300, 300,  80, 300, 190)
x3 <- c( 10,  10, 120, 120,  10,  10,  65,  65,  54, 120, 120,  10, 120)
rate <- c(8.55,  3.79, 4.82, 0.02,  2.75, 14.39, 2.54,
          4.35, 13.00, 8.50, 0.05, 11.32,  3.13)
fun <- function(b)
        (b[1]*x2 - x3/b[5])/(1 + b[2]*x1 + b[3]*x2 + b[4]*x3) - rate
lsqnonlin(fun, rep(1, 5))
# $x    [1.25258502 0.06277577 0.04004772 0.11241472 1.19137819]
# $ssq  0.298901

##  Example for lsqnonneg()
C <- matrix(c(0.0372, 0.2868,
              0.6861, 0.7071,
              0.6233, 0.6245,
              0.6344, 0.6170), nrow = 4, ncol = 2, byrow = TRUE)
d <- c(0.8587, 0.1781, 0.0747, 0.8405)

sol <- lsqnonneg(C, d)
cbind(qr.solve(C, d), sol$x)
# -2.563884  5.515869e-08
#  3.111911  6.929003e-01

##  Example for lsqcurvefit()
#   Lanczos1 data (artificial data)
#   f(x) = 0.0951*exp(-x) + 0.8607*exp(-3*x) + 1.5576*exp(-5*x)
x <- linspace(0, 1.15, 24)
y <- c(2.51340000, 2.04433337, 1.66840444, 1.36641802, 1.12323249, 0.92688972,
       0.76793386, 0.63887755, 0.53378353, 0.44793636, 0.37758479, 0.31973932,
       0.27201308, 0.23249655, 0.19965895, 0.17227041, 0.14934057, 0.13007002,
       0.11381193, 0.10004156, 0.08833209, 0.07833544, 0.06976694, 0.06239313)

p0 <- c(1.2, 0.3, 5.6, 5.5, 6.5, 7.6)
fp <- function(p, x) p[1]*exp(-p[2]*x) + p[3]*exp(-p[4]*x) + p[5]*exp(-p[6]*x)
lsqcurvefit(fp, p0, x, y)

##  Example for lsqsep()
f <- function(x) 0.5 + x^-0.5 + exp(-0.5*x)
set.seed(8237); n <- 15
x <- sort(0.5 + 9*runif(n))
y <- f(x)                       #y <- f(x) + 0.01*rnorm(n)

m <- 2
f1 <- function(b, x) x^b
f2 <- function(b, x) exp(b*x)
flist <- list(f1, f2)
start <- c(-0.25, -0.75)

sol <- lsqsep(flist, start, x, y, const = TRUE)
a0 <- sol$a0; a <- sol$a; b <- sol$b
fsol <- function(x) a0 + a[1]*f1(b[1], x) + a[2]*f2(b[2], x)

\dontrun{
    ezplot(f, 0.5, 9.5, col = "gray")
    points(x, y, col = "blue")
    xs <- linspace(0.5, 9.5, 51)
    ys <- fsol(xs)
    lines(xs, ys, col = "red")
}
}
\keyword{ fitting }
back to top