Using different or custom root solvers in geex

Choice of rootFUN

By default, geex uses the rootSolve::multiroot function for finding roots of a set of estimating equations when compute_roots = TRUE in m_estimate(). However, a user can choose a different root find find algorithm via the root_control argument.

For example, consider the following estFUN which is Huber’s estimator for the center of symmetric distributions [@stefanski2002; example 6]. This example was chosen because it has a single root, so that the stats::uniroot function can be used to find the roots.

myefun <- function(data, k = 1.5){
  function(theta){
    x <- data$Y1 - theta[1]
    if(abs(x) <= k) x else sign(x) * k
  }
}

Internally, estFUN is used to build $G_m = \sum_{i = 1}^m \psi(O_i, \theta)$ or in psuedo-code f = sum(inner_estFUN(theta)). f is passed to the root finding function along with options in the root_control arguments. For example, multiroot requires f and start (starting values for the algorithm:

multiroot_results <- m_estimate(
  estFUN = myefun, 
  data  = geexex,
  root_control = setup_root_control(start = 3))

The stats::uniroot function, however, requires the arguments f and interval (or lower and upper)

uniroot_results <- m_estimate(
  estFUN = myefun, 
  data  = geexex,
  root_control = setup_root_control(stats::uniroot, interval = c(0, 10)))

Comparing results:

roots(multiroot_results) - roots(uniroot_results)
## [1] 6.244845e-07

They are basically the same, but this may not be true depending f and the options given to the root finder.

Building a custom root finder

All that is necessary for rootFUN is a function where:

  • the first argument is the function whose roots are sought
  • the output is a named list where the root estimates can be identified by the rootFUN_object argument in m_estimate. For example, both uniroot and multiroot return a list where the root estimates are in the item named “roots”. The default is rootFUN_object = 'roots', so this option works for both uniroot and multiroot.