update(fit) on sdmSimple uses wrong link function for c

[venpopov]

For computational efficiency, I had specified c to have an identity link in the brms custom family, but coded it manually in the likelihood function on the log scale (relative to the model equations).

To not confuse users, I change the link to log in the postprocessing method postprocess_brm.sdmSimple so that it displays as a log-link in the summary(fit).

However, I just realized that then if you update the model using

fit_new <- update(fit)

brms now will think it needs to use a log link and essentially apply an extra exp transformation.

Two solutions:

• give up the manual coding of the log link at the cost of some computation time (i have to check how much, maybe it's not a huge loss)
• create an update.bmmmodel method, and assign class(fit) to be bmmmodel, in which we take care of that

Related issues

We also had the idea to hide irrelevant terms of the formula from the summary. This will present the same problem, as then the formula would be changed for the update function.

[GidonFrischkorn]

For the SDM parameterization, an alternative could be to explicitly state that the parameters are logC and kappa. Then you could add in the information on the model that C is parameterized on the log scale for computational efficiency, but technically has an identity link. Maybe we can ask Klaus if he thinks that would be confusing.

With respect to cleaning up the summary, maybe we could also have a look into the summary.brmsfit method and can see if we can adapt this method to a summary.bmmfit that prints only the elements of the function provided by the user.

[venpopov]

A few corrections to what we discussed in person.

If we have an update method like this (with more in it):

#' @export
update.bmmfit <- function(fit) {
  message('Yay!')
  fit = NextMethod(fit)
  fit$testing = 'abc'
  message('Also yay!')
  return(fit)
}

And do

fit <- fit_model(bmf(c ~ 0 + cond, kappa ~ 0 + cond),
                 data = dat,
                 model = sdmSimple('y'),
                 parallel = T,
                 chains = 4,
                 iter=500,
                 file = 'local/fit_sdm_m',
                 backend='cmdstanr')

class(fit) <- c('bmmfit','brmsfit')
fit <- update(fit)

It actually does return the fit with the postprocessing done after NextMethod, in this case having fit$testing = 'abc'. When we looked at it in the office, I just ran update(fit) but did not assign the result to the fit variable, so we were looking at the old fit object... so in principle we can make update wrapper, contrary to what we thought. Let's leave this issue open as a todo

[venpopov]

I will actually leave this in the code now, both as a reminder of how it will work, and to prevent users from trying to use update:

#' @export
update.bmmfit <- function(fit, ...) {
  stop("The update method for bmmfit is not yet implemented, but is planned for a future release.", call. = FALSE)
  # do any necessary preprocessing to accomondate bmm changes into brms
  # ....

  fit = NextMethod(fit)  # pass back to brms::update.brmsfit

  # do any necessary postprocessing to convert back to bmmfit
  # ....
  class(fit) <- c('bmmfit', class(fit))  # reassing class
  return(fit)
}