daarem {daarem} | R Documentation |
An ‘off-the-shelf’ acceleration scheme for accelerating the convergence of any smooth, monotone, slowly-converging fixed-point iteration. It can be used to accelerate the convergence of a wide variety of montone iterations including, for example, expectation-maximization (EM) algorithms and majorization-minimization (MM) algorithms.
daarem(par, fixptfn, objfn, ..., control=list())
par |
A vector of starting values of the parameters. |
fixptfn |
A vector function, |
objfn |
This is a scalar function, |
control |
A list of control parameters specifying any changes to default values of algorithm control parameters. Full names of control list elements must be specified, otherwise, user-specifications are ignored. See *Details*. |
... |
Arguments passed to |
Default values of control
are:
maxiter=2000
,
order=10
,
tol=1e-08
,
mon.tol=0.01
,
cycl.mon.tol=0.0
,
alpha=1.2
,
kappa=25
,
resid.tol=0.95
,
convtype="param"
maxiter
An integer denoting the maximum limit on the number
of evaluations of fixptfn
, G
. Default value is 2000.
order
An integer 1
denoting the order of the DAAREM acceleration scheme.
tol
A small, positive scalar that determines when iterations
should be terminated. When convtype
is set to "param", iteration is terminated when
||x_k - G(x_k)|| < tol
.
Default is 1.e-08
.
mon.tol
A nonnegative scalar that determines whether the montonicity condition
is violated. The monotonicity condition is violated whenver L(x[k+1]) < L(x[k]) - mon.tol
.
Such violations determine how much damping is to be applied on subsequent steps of the algorithm. Default
value of mon.tol is 1.e-02
.
cycl.mon.tol
A nonegative scalar that determines whether a montonicity condition is violated
after the end of the cycle. This cycle-level monotonicity condition is violated whenver
L(x[end cycle]) < L(x[start cycle]) - cycl.mon.tol
.
Here, x[start cycle]
refers to the value of x
at the beginning of the current cycle
while x[end cycle]
refers to the value of x
at the end of the current cycle.
Such violations also determine how much damping is to be applied on subsequent steps of the algorithm.
kappa
A nonnegative parameter which determines the “half-life” of relative damping and how quickly relative
damping tends to one. In the absence of monotonicity
violations, the relative damping factor is <= 1/2
for the first kappa
iterations, and it is
then greater than 1/2
for all subsequent iterations. The relative damping factor is the ratio between
the norm of the unconstrained coefficients in Anderson acceleration and the norm of the damped coefficients.
In the absence of any monotonicity violations, the relative damping factor in iteration k
is
1/(1 + \alpha^(\kappa - k))
.
alpha
A parameter > 1
that determines the initial relative damping factor and how quickly
the relative damping factor tends to one. The initial relative damping factor is 1/(1 + \alpha^\kappa)
.
In the absence of any monotonicity violations, the relative damping factor in iteration k
is
1/(1 + \alpha^(\kappa - k))
.
resid.tol
A nonnegative scalar < 1
that determines whether a residual change condition is violated.
The residual change condition is violated whenever ||x_k+1 - G(x_k+1)|| > ||x_k - G(x_k)|| (1 + resid.tol^k)
. Default value of resid.tol is 0.95
.
convtype
This can equal either "param" or "objfn". When set to "param", convergence is determined by the criterion: ||x_k - G(x_k)|| \leq tol
. When set to "objfn", convergence is determined by the objective function-based criterion: | L(x[k+1]) - L(x[k])| < tol
.
intermed
A logical variable indicating whether or not the intermediate results of iterations should be returned. If set to
TRUE
, the function will return a matrix where each row corresponds to parameters at each iteration, along with the corresponding value of the objective function in the first column. This option is inactive when objfn is not specified. Default is FALSE
.
A list with the following components:
par |
Parameter, |
value.objfn |
The value of the objective function |
fpevals |
Number of times the fixed-point function |
objfevals |
Number of times the objective function |
convergence |
An integer code indicating type of convergence. |
objfn.track |
A vector containing the value of the objective function at each iteration. |
p.intermed |
A matrix where each row corresponds to parameters at each iteration,
along with the corresponding value of the objective function (in the first column).
This object is returned only when the control parameter |
Nicholas Henderson and Ravi Varadhan
Henderson, N.C. and Varadhan, R. (2019) Damped Anderson acceleration with restarts and monotonicity control for accelerating EM and EM-like algorithms, Journal of Computational and Graphical Statistics, Vol. 28(4), 834-846. doi: 10.1080/10618600.2019.1594835
n <- 2000
npars <- 25
true.beta <- .5*rt(npars, df=2) + 2
XX <- matrix(rnorm(n*npars), nrow=n, ncol=npars)
yy <- ProbitSimulate(true.beta, XX)
max.iter <- 1000
beta.init <- rep(0.0, npars)
# Estimating Probit model with DAAREM acceleration
aa.probit <- daarem(par=beta.init, fixptfn = ProbitUpdate, objfn = ProbitLogLik,
X=XX, y=yy, control=list(maxiter=max.iter))
plot(aa.probit$objfn, type="b", xlab="Iterations", ylab="log-likelihood")
# Compare with estimating Probit model using the EM algorithm
max.iter <- 25000 # need more iterations for EM convergence
beta.init <- rep(0.0, npars)
em.probit <- fpiter(par=beta.init, fixptfn = ProbitUpdate, objfn = ProbitLogLik,
X=XX, y=yy, control=list(maxiter=max.iter))
c(aa.probit$fpevals, em.probit$fpevals)
c(aa.probit$value, em.probit$value)
# Accelerating using SQUAREM if the SQUAREM package is loaded
# library(SQUAREM)
# max.iter <- 5000
# sq.probit <- squarem(par=beta.init, fixptfn=ProbitUpdate, objfn=ProbitLogLik,
# X=XX, y=yy, control=list(maxiter=max.iter))
# print( c(aa.probit$fpevals, em.probit$fpevals, sq.probit$fpevals) )
# print( c(aa.probit$value, em.probit$value, sq.probit$value) )
# print( c(aa.probit$objfeval, em.probit$objfeval, sq.probit$objfeval) )