Please also provide the following information in addition to your question:

- Operating System: Linux Ubuntu
- brms Version: 2.10.0

When developing a new product blend, I run a response surface design of experiment. The results are then analyzed using linear regression. Often I vary e.g four factors at the time, and the design will support main, interaction, and squared effects = Intercept + 4 main + 4 squared + 6 interaction terms.

Features are centered and scaled.

In a frequentist setting, I would use leave one out cross validation and remove terms that improve predictive performance (Q2 = 1 - (prediction residual sum of squares / total sum of squares)).

Terms are removed observing a hierarchy, meaning I cannot leave out a main effect but still leave in an interaction term or a squared effect.

I would like to build an automatic pipeline to analyze these kinds of experiments, and I would like to use a Bayesian approach, preferable using brms. Now Iâ€™m wondering how to select the best subset in a Bayesian setting using brms.

- shrinkage: should I opt for shrinkage using either a horsehoe prior on the effects or use a multilevel model, assuming that effects are normally distributed? Using this route, I do not need to worry about feature hierarchy! However, will there be enough predictors (e.g. 8+) in my setting for this to make sense?
- fit several models, representing all valid formulas and find the best performing model according LOO?
- something entirely differentâ€¦

If 2) is the best option, do you have any ideas how to automatic generate all valid formulas according to the hierarchy criteria outlined above?

Any insights, ideas, and thoughts are most welcome.

Thank you for being e great community.

/Jannik