9
$\begingroup$

I have a dataset with three variables, where all variables are quantitatives. Let call it $y$, $x_1$ and $x_2$. I'm fitting a regression model in a Bayesian perspective via MCMC with rjags

I done a exploratory analysis and the scatterplot of $y\times x_2$ suggest that a quadratic term should be used. Then I fitted two models

(1) $y=\beta_0+\beta_1*x_1+\beta_2*x_2$

(2) $y=\beta_0+\beta_1*x1+\beta_2*x_2+\beta_3*x_1x_2+\beta_4*x_1^2+\beta_5*x_2^2$

In model 1 the effect size of each parameter is not small and the 95% credible interval not contains the value $0$.

In model 2 the effect size of parameters $\beta_3$ and $\beta_4$ are small and each of credible intervals for all parameters contains $0$.

The fact that a credible interval contains $0$ is enough to say that the parameter is not significant?

Then I adjusted the following model

(3)$y=\beta_0+\beta_1*x_1+\beta_2*x_2+\beta_3*x^2_2$

The effect size of each parameter is not small, but with exception of $\beta_1$ all credible intervals contains $0$.

Which is the right way to do variable selection in Bayesian statistics?

EDIT: I can use Lasso in any regression model, like Beta model? I'm using a model with variable dispersion where $$log(\sigma)=-\pmb{\delta}X$$ where $\pmb{\delta}$ is a vector. I should use Laplace prior in $\pmb{\delta}$ too?

EDIT2: I fitted two models, one with Gaussian priori for $\beta_j$, $\delta_j$ and one with Laplace(double-exponential).

The estimatives for the Gaussian model are

            Mean      SD  Naive SE Time-series SE
B[1]     -1.17767 0.07112 0.0007497      0.0007498
B[2]     -0.15624 0.03916 0.0004128      0.0004249
B[3]      0.15600 0.05500 0.0005797      0.0005889
B[4]      0.07682 0.04720 0.0004975      0.0005209
delta[1] -3.42286 0.32934 0.0034715      0.0034712
delta[2]  0.06329 0.27480 0.0028966      0.0028969
delta[3]  1.06856 0.34547 0.0036416      0.0036202
delta[4] -0.32392 0.26944 0.0028401      0.0028138

The estimatives for the Lasso model are

              Mean      SD  Naive SE Time-series SE
B[1]     -1.143644 0.07040 0.0007421      0.0007422
B[2]     -0.160541 0.05341 0.0005630      0.0005631
B[3]      0.137026 0.05642 0.0005947      0.0005897
B[4]      0.046538 0.04770 0.0005028      0.0005134
delta[1] -3.569151 0.27840 0.0029346      0.0029575
delta[2] -0.004544 0.15920 0.0016781      0.0016786
delta[3]  0.411220 0.33422 0.0035230      0.0035629
delta[4] -0.034870 0.16225 0.0017103      0.0017103
lambda    7.269359 5.45714 0.0575233      0.0592808

The estimatives for $\delta_2$ and $\delta_4$ reduced a lot in Lasso model, it means that I should remove this variables from the model?

EDIT3: The model with double exponential prior (Lasso) gaves me bigger Deviance, BIC and DIC values than the model with Gaussian priors and I even get a smaller values after removing the dispersion coefficient $\delta_2$ in the Gaussian model.

$\endgroup$
  • 2
    $\begingroup$ Section 18.4 of DBDA2E* is about variable selection in multiple regression. With great caution you can put in inclusion indicators with each coefficient and look at the posterior probability of inclusion. When interpreting posterior distributions of parameters, the 95%HDI including zero does not indicate equivalence to zero. *DBDA2E = Doing Bayesian Data Analysis 2nd Edition. $\endgroup$ – John K. Kruschke May 17 '17 at 22:55
  • 2
    $\begingroup$ The natural way to compare models in a Bayesian framework is through marginal likelihoods, not credible intervals. An alternative connected with model averaging is to use a mixture representation and infer from the weights of each model/component which model is favoured by the data. $\endgroup$ – Xi'an May 19 '17 at 7:43
  • $\begingroup$ @Xi'an but comparing two or more models through marginal likelihoods wouldn't be the same as using Bayes factors, if all models were given the same prior probability? $\endgroup$ – DeltaIV May 19 '17 at 17:48
  • $\begingroup$ Dear Prof Kruschke, I have a doubt regarding the computation of crdible intervals. What I have understood that there can be many credible intervals depending upon plausibility of posterior based upon different priors. But here how to decide which priors are most plausible which in turn gives most plausible posterior? Another question is, I am using Variational Inference (VI) for computing the posteriors and eventually compute the lower bound of model evidence. How to compute the credible interval for the posteriors in case of VI? Moreover how to proceed for the Bayes Factor in case of VI? $\endgroup$ – Sandipan Karmakar Sep 15 '17 at 9:28
9
+100
$\begingroup$

It is well known that building a model based on what is significant (or some other criterion such as AIC, whether a credible interval contains 0 etc.) is pretty problematic, particularly if you then do inference as if you had not done model building. Doing a Bayesian analysis does not change that (see also https://stats.stackexchange.com/a/201931/86652). I.e. you should not do variable selection, but rather model averaging (or something that could get you some zero coefficients, but reflects the whole modelling process, such as LASSO or elastic net).

Bayesian model choice is more typically framed as Bayesian model averaging. You have different models, each with a different prior probability. If the posterior model probability for a model becomes low enough, you are essentially entirely discarding the model. For equal prior weights for each model and flat priors, model averaging with weights proportional to $\exp(-\text{BIC}/2)$ for each model approximates this.

You can alternatively express the model averaging as a prior that is a mixture between a point mass (the weight of the point mass is the prior probability of the effect being exactly zero = the effect is not in the model) and a continuous distribution (e.g. spike-and-slab priors). MCMC sampling can be quite difficult for such a prior.

Carvalho et al. motivate the horseshoe shrinkage prior by suggesting that it works like a continuous approximation to a spike-and-slab prior. It is also a case of embedding the problem in a hierarchical model, where to some extent the size and presence of effects on some variables relax the required evidence for others a bit (through the global shrinkage parameter, this is a bit like false-discovery rate control) and on the other hand allow individual effects to stand on their own if the evidence is clear enough. There is a convenient implementation of it available from the brms R package that builds on Stan/rstan. There are a number of further similar priors such as the horseshoe+ prior and the whole topic is an area of ongoing research.

$\endgroup$
  • $\begingroup$ Bayesian Lasso is like this stats.stackexchange.com/questions/28609/… ? I'm a model with dispersion variable I should use the double exponentail prior for those parameters too? $\endgroup$ – user72621 May 19 '17 at 13:58
2
$\begingroup$

There are a number of formal methods for Bayesian variable selection. A slightly outdated review of Bayesian variable selection methods is presented in:

A review of Bayesian variable selection methods: what, how and which

A more recent review, which also includes a comparison of different methods and the performance of R packages where they are implemented is:

Methods and Tools for Bayesian Variable Selection and Model Averaging in Univariate Linear Regression

This reference is particularly useful in that it points you to specific R packages where you just need to plug in the response and the covariate values (and in some cases the hyperparameter values) in order to run the variable selection.

Another, quick and dirty and non-recommended, way of conducting "Bayesian" variable selection is to use stepwise selection (forward, backward, both) using BIC and the R command stepAIC(), which can be tweaked to perform selection in terms of BIC.

https://stat.ethz.ch/R-manual/R-devel/library/MASS/html/stepAIC.html

Another quick and dirty way of testing $\beta_4=0$ is by using the Savage-Dickey density ratio and the posterior simulation you already got:

https://arxiv.org/pdf/0910.1452.pdf

$\endgroup$
  • $\begingroup$ I think the question is asking why the three parameters in model 3 all have credible regions containing 0 and not whether or not $\beta_4$ is 0. $\endgroup$ – Michael R. Chernick May 19 '17 at 15:44
  • $\begingroup$ @MichaelChernick Then why is the OP asking "In this case is reasonable say that $\beta_4\neq 0$"? and "Which is the right way to do variable selection in Bayesian statistics"? $\endgroup$ – CTHULHU May 19 '17 at 15:48
  • $\begingroup$ I missed that part of the question but I don't think it was the main issue. $\endgroup$ – Michael R. Chernick May 19 '17 at 15:50
  • 1
    $\begingroup$ @MichaelChernick Well, I guess the OP has the last word here ... $\endgroup$ – CTHULHU May 19 '17 at 15:51
1
$\begingroup$

The whole idea of Bayesian statistics is different from a frequentist approach. In this way I think to use the terms of significance is not accurate. I guess it is up to the reader to decide if the results (distribution) you get from your model for your $β$ 's are for him reliable or trustful. It always depends on the distribution itself. How skewed and wide is it and how much of the area is below zero?

You can also find a nice lecture about the topic here at 41:55 :

https://vimeo.com/14553953

$\endgroup$
  • $\begingroup$ I added a example of histogram of one variable which credible interval contains $0$ could you take a look? $\endgroup$ – user72621 May 19 '17 at 15:17
  • 1
    $\begingroup$ Back from the weekend. Where can we find the histogram? $\endgroup$ – burton030 May 22 '17 at 7:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy