Predict.brmsfit fails with a model that has a monotonic predictor

Interfaces brms
1

I am running the following model that has one numeric predictor, x1, and one ordered factor monotonic predictor, x2:

f <- as.formula(y ~ x1 + mo(x2))
model <- brm(f, data=mydata,iter=9000,control=list(adapt_delta=0.9))

However, when I try to enter it into posterior_predict:

myprediction <- as.data.frame(brms::posterior_predict(model,newdata=mydata))

I receive the following error:

Error in simplex[, X + 1] : subscript out of bounds

I am able to run posterior_predict when I remove mo() from the formula or when I convert the ordered factor x2 to numeric --which makes me believe the monotonic predictor is the problem. Am I using mo() the wrong way, or are these two functions simply incompatible? Thanks for your advice.

2

Looks like a bug. Let me check.

3

Can you please provide a minimal reproducible example, ideally with simulated data? I don’t have your data so I cannot run your code.

4

Hi @paul.buerkner I believe I have figured out how to fix the problem, although I’m still not entirely sure why the issue occurs. I’ve posted a reproducible example below with explanation in case other users have the same problem. To make the sample data similar to my data, I added a few random NA values into the factor column.

Normally I don’t need to manually exclude NA values from my data when running brms models because the model automatically excludes these when run; and indeed when I convert my ordered factor variable into numeric, posterior_predict works fine with NA values included in the input dataset. However, there seems to be a special case where posterior_predict fails to evaluate predictions for a dataset that has NA values in a factor using a model where NA values had been automatically excluded as a process of the brm function.

For example, the following code returns my original error:

# R version 4.0.2
library(brms) # version brms_2.14.4
library(ggplot2) # version ggplot2_3.3.2
library(dplyr) # version dplyr_1.0.2
mydata <- mpg
mydata$cyl <- factor(mydata$cyl,ordered=TRUE) 
mydata <- mydata %>% mutate_if(is.numeric,scale)
mydata$cyl[10:15] <- NA

f <- as.formula(hwy ~ displ + mo(cyl))
model <- brm(f, data=mydata,iter=9000,control=list(adapt_delta=0.9))
myprediction <- as.data.frame(brms::posterior_predict(model,newdata=mydata))

But when I add mydata <- mydata[complete.cases(mydata),] before running the model the error disappears.