# Unable to fit nonlinear robit model to new data

**URL:** <https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284>\
**Category:** brms\
**Created:** [August 6, 2020, 5:39pm UTC](https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284 "2020-08-06T17:39:41Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![josh\_fiechter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/josh_fiechter/32/6031_2.png) [@josh\_fiechter](https://discourse.mc-stan.org/u/josh_fiechter)\
**Post date:** [August 6, 2020, 5:39pm UTC](https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284/1 "2020-08-06T17:39:41Z")

</div>

I fit a nonlinear robit model in `brms` that was set up as follows:

```
stan_inv_robit <- "
  real inv_robit(real y, real nu) {
    return(student_t_cdf(y, nu, 0, (nu - 2) / nu));
  }
"
stanvar_inv_robit <- stanvar(scode = stan_inv_robit, block = "functions")

bf_1 <- bf(Accuracy_rev ~ inv_robit(eta, nu),
           nu ~ 1,
           b + m ~ 1,
           bl + tau ~ 1 + (1|Subset + User + Token),
           nlf(eta ~ (tau-((bl + repetition)^.1)*(model_time_z^-(b + m * stability_z)))/.1),
           nl = T,
           family = bernoulli('identity'))

prior_1 <-
  prior(gamma(2, 0.1), nlpar = "nu", lb = 2) +
  prior(normal(0, .5), nlpar = 'b') +
  prior(normal(0, .5), nlpar = 'm') +
  prior(normal(0, .5), nlpar = 'bl') +
  prior(normal(0, .5), nlpar = 'tau')

fit_1 <- brm(bf_1,
             prior = prior_1,
             stanvars = stanvar_inv_robit,
             data = filter(dat_prepped, `Train or Test` == 1),
             chains = 4)

```

I can obtain model fits for the data set that produced the parameter estimates, but when I try to get fitted values for a new set of data (i.e., `fitted(fit_1, newdata = test, allow_new_levels = T)` I get the following error:

`Error: Error in (function (y, nu, pstream__ = <pointer: 0x000000000211ab40>) : Exception: student_t_cdf: Random variable is nan, but must not be nan! (in 'unknown file name' at line 5)`

At first I thought my test set might have been returning `NA` values of `eta` that are passed on to the `inv_robit` function, but I was able to calculate non-`NA` values for every observation using the estimated parameter values. So I’m not sure what the issue is. Thanks for any advice!

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [August 6, 2020, 6:37pm UTC](https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284/2 "2020-08-06T18:37:31Z")

</div>

Can you show all the post processing code please? It seems you only showed parts of the code that is relevant.

---

<div class="post-metadata">

**Author:** ![josh\_fiechter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/josh_fiechter/32/6031_2.png) [@josh\_fiechter](https://discourse.mc-stan.org/u/josh_fiechter)\
**Post date:** [August 6, 2020, 8:29pm UTC](https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284/3 "2020-08-06T20:29:32Z")

</div>

Sure – let me know if there’s anything else that might be helpful.

```
# Expose inv_robit function
expose_functions(fit_1, vectorize = T)

# Get training observations from main data set
train <- na.omit(filter(dat_prepped, `Train or Test` == 1))
train$fitted <- fitted(fit_1)[,1] # Runs fine

# Get test observations from main data set
test <- na.omit(filter(dat_prepped, `Train or Test` == 0))
test$fitted <- fitted(fit_1, newdata = test, allow_new_levels = T)[,1] # Generates an error

# Get values of eta from model parameters and see if any are NA
test <- mutate(test, b = fitted(fit_1, newdata = test, allow_new_levels = T, nlpar = 'b')[,1],
               m = fitted(fit_1, newdata = test, allow_new_levels = T, nlpar = 'm')[,1],
               bl = fitted(fit_1, newdata = test, allow_new_levels = T, nlpar = 'bl')[,1],
               tau = fitted(fit_1, newdata = test, allow_new_levels = T, nlpar = 'tau')[,1],
               eta = (tau-((bl + repetition)^.1)*(model_time_z^-(b + m * stability_z)))/.1)

length(which(is.na(test$eta))) # Returns 0

```

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [August 6, 2020, 9:36pm UTC](https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284/4 "2020-08-06T21:36:19Z")

</div>

Hmm unfortunately nothing in that code to identify the problem right away. Could you make a minimal reproducible example so that I can check myself?

---

<div class="post-metadata">

**Author:** ![josh\_fiechter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/josh_fiechter/32/6031_2.png) [@josh\_fiechter](https://discourse.mc-stan.org/u/josh_fiechter)\
**Post date:** [August 9, 2020, 7:22pm UTC](https://discourse.mc-stan.org/t/unable-to-fit-nonlinear-robit-model-to-new-data/17284/5 "2020-08-09T19:22:20Z")

</div>

While in the process of making up an example, I realized that I was getting `NaN` values for `eta` if I used `summary = F` in the `fitted` function. And that seems to be because my priors weren’t reflecting the fact that I’m using the identity link. So this issue seems to have been completely my fault. I’ll re-run with better priors and let you know if anything pops up.
