pbs-assess / sdmTMB

:earth_americas: An R package for spatial and spatiotemporal GLMMs with TMB
https://pbs-assess.github.io/sdmTMB/
186 stars 26 forks source link

Add argument to _mix() families for fixed proportion? #318

Open seananderson opened 7 months ago

seananderson commented 7 months ago

We should do some testing but I'm not sure you can actually estimate the proportion very often. The Laplace approximation may not work well here.

It can currently be fixed by doing something like this

fit <- sdmTMB(
  ...,
  family = delta_lognormal_mix(), 
  control = sdmTMBcontrol(
    start = list(logit_p_mix = qlogis(0.01)), # 1 % extreme distribution
    map = list(logit_p_mix = factor(NA))
  )
)

but that's pretty clunky if you almost always need to do it.

We may want to disable estimation of this proportion entirely.

ericward-noaa commented 6 months ago

I'm not sure what situations this is going to be well estimated. I had done some initial tests with examples like this data rich situation where it can recover the log_ratio_mix parameter -- but haven't explored less data / larger ranges of the mixture parameter / fraction of zeros / etc

set.seed(123)
  range <- 1
  x <- stats::runif(5000, -1, 1)
  y <- stats::runif(5000, -1, 1)
  loc <- data.frame(x = x, y = y)
  spde <- make_mesh(loc, c("x", "y"), n_knots = 70, type = "kmeans")
  sigma_O <- 0.3
  sigma_E <- 0
  phi <- 0.2
  s <- sdmTMB_simulate(~ 1, loc, mesh = spde, family = lognormal(),
                       B = 1,
                       phi = phi, range = range, sigma_O = sigma_O, seed = 1
  )

  zeros <- sample(1:nrow(s), size=2000, replace=F)
  second_mix <- sample(1:nrow(s),size=500, replace=F)

  s$observed[second_mix] <- exp(log(s$observed[second_mix]) + 2)
  s$observed[zeros] <- 0
  mlog <- sdmTMB(data = s, formula = observed ~ 1, mesh = spde,
                 family = delta_lognormal_mix())

  expect_equal(log(exp(mlog$sd_report$par.fixed[["log_ratio_mix"]])+1), 2.0, tolerance = 0.01)