A message from our CEO about the future of Stack Overflow and Stack Exchange. Read now.

Hot answers tagged

19

The beta regression approach is to reparameterize in terms of $\mu$ and $\phi$. Where $\mu$ will be the equivalent to y_hat that you predict. In this parameterization you will have $\alpha=\mu\times\phi$ and $\beta=(1-\mu) \times \phi$. Then you can model $\mu$ as the logit of the linear combination. $\phi$ can either have its own prior (must be greater ...


18

Greg Snow gave a great answer. For completeness, here is the equivalent in Stan syntax. Although Stan has a beta distribution that you could use, it is faster to work out the logarithm of the beta density yourself because the constants log(y) and log(1-y) can be calculated once at the outset (rather than every time that y ~ beta(alpha,beta) would be called). ...


18

From (an earlier version of) the Stan reference manual: Not specifying a prior is equivalent to specifying a uniform prior. A uniform prior is only proper if the parameter is bounded[...] Improper priors are also allowed in Stan programs; they arise from unconstrained parameters without sampling statements. In some cases, an improper prior may ...


18

Firstly, feel free to ask questions like this on our users' list (http://mc-stan.org/mailing-lists.html) where we discuss not only issues related to Stan implementations/optimizations/etc but also practical statistical and modeling questions. As to your question, it's absolutely a fine approach. There are many ways to justify it more rigorously (for ...


16

Quick thoughts: 1) The key issue is what applied question you are trying to answer for your audience, because that determines what information you want from your statistical analysis. In this case, it seems to me that you want to estimate the magnitude of differences between groups (or perhaps the magnitude of ratios of the groups if that is the measure ...


12

Because all of the parameters of this distribution are known, and we merely want to draw samples from this distribution, coding the model in rstan is straightforward. Note that this is, by far, one of the least efficient paths to sampling from this particular model, in terms of time that I spent coding it (15 minutes). The author of the original post is ...


12

This is of course a diverse set of people with a range of opinions getting together and writing a wiki. I summarize I know/understand with some commentary: Choosing your prior based on computational convenience is an insufficient rationale. E.g. using a Beta(1/2, 1/2) solely because it allows conjugate updating is not a good idea. Of course, once you ...


9

I would distinguish between the prior distribution and the parameters for the prior distribution. When I started with Stan, I would set the parameters to the prior distributions just as some values. So in the model step, I would have something like model { mu ~ normal(0, 1) y ~ normal(mu, s) } for a normal prior on the mean coefficient for ...


9

Hamiltonian Monte Carlo is often quoted as being rotation invariant, but what does that actually mean? In theory Hamiltonian Monte Carlo (HMC) is independent of the coordinates chosen to represent your distribution, which actually means it is both rotation and scale invariant. The problem is that, in practice, we can't run HMC exactly and must instead use ...


9

Fitting random slopes with the population-level slope fixed to zero is not out of the question - it's not mathematically or statistically ill-posed - but it's a rather weird model that would require some extra justification. Why would you expect that the average slope across cities would be exactly zero (which is what is implied by the model that omits the ...


8

From the Stan reference v1.0.2 (pg 6, footnote 1) If no prior were specified in the model block, the constraints on theta ensure it falls between 0 and 1, providing theta an implicit uniform prior. For parameters with no prior specified and unbounded support, the result is an improper prior. Stan accepts improper priors, but posteriors must be proper in ...


8

Keep in mind that in Bayesian inference, you are looking at the entire range of values in the posterior, rather than point-estimate summaries. After convergence, MCMC samples are samples drawn from that posterior density. Using just the latest sample is the same as picking a random point from the posterior and using that as your sample estimate. For an ...


8

As my previous answer was deleted, here is a more explicit one, with an example using a simple sampling from the prior: library(rstan) model = " parameters { real p; } model { p ~ normal(1,3); } " fit = stan(model_code = model, pars = c('p'), control=list(adapt_delta=0.99, max_treedepth=10), iter = 5000, chains = 1, ...


8

If you have a disjoint set of possible events (say three of them) that exhaust all possiblites, the probabilities of each of these events must sum to one: $$p_1 + p_2 + p_3 = 1$$ Being probabilities, each of these is also bounded between zero and one: $$ 0 \leq p_i \leq 1 $$ If consider the three probabilities as a point in euclidean space $(p_1, p_2, ...


8

They do not provide any scientific/mathematical justification for doing so. Most of the developers do not work on this kind of priors, and they prefer to use more pragmatic/heuristic priors, such as normal priors with large variances (which may be informative in some cases). However, it is a bit strange that they are happy to use PC priors, which are based ...


7

You can define a proper or improper prior in the Stan language using the increment_log_prob() function, which will add its input to the accumulated log-posterior value that is used in the Metropolis step to decide whether to accept or reject a proposal for the parameters. In your example, the model block would need to include the new line ...


7

As mentioned by previous answers, Stan, JAGS, and WinBUGS require that priors be specified as mathematical functions. If you've already got an MCMC-represented posterior from a previous analysis, and you want to use that MCMC posterior as a prior for subsequent data, you must approximate the MCMC posterior in a mathematical form. Unless you have a simple ...


7

Likelihood For a mixture of two Gaussians, the likelihood can be written as: $$ y_i \sim \pi N(y_i|\alpha_0 + x_i\beta, \sigma_0) + (1-\pi) N(y_i|\alpha_1 + x_i\beta, \sigma_1) $$ where $\pi \in [0, 1]$. This is fine, but having two components in the likelihood makes sampling more difficult. A trick when dealing with mixture models is to augment the ...


7

A divergent transition in Stan tells you that the region of the posterior distribution around that divergent transition is geometrically difficult to explore. For example here is a quote from the manual: The primary cause of divergent transitions in Euclidean HMC (other than bugs in the code) is highly varying posterior curvature, for which small step ...


6

$\hat{R}$ and "potential scale reduction factor" refer to the same thing. See Chapter 6 of the Handbook of Markov Chain Monte Carlo, "Inference from Simulations and Monitoring Convergence" by Andrew Gelman and Kenneth Shirley. In Stan, the number reported is actually split $\hat{R}$; the calculation of $\hat{R}$ is computed with each of the chains split in ...


6

The primary reason that your code does not yield the expected answer is that you are using the multi_normal_prec likelihood rather than the multi_normal likelihood. The former expects a precision matrix (the inverse of a covariance matrix) as its second argument, while the latter expects a covariance matrix. For what it is worth, you should be able to ...


6

The main question was the result of a bug in rstanarm that has since been fixed on GitHub. However, in general, we do not recommend rstanarm models that exclude the intercept. A better alternative is to place a tight prior with mean zero on the intercept. In this case, including the intercept yields a better fit and similar results between lme4::lmer and ...


6

I'll leave this as an "answer" as I don't have enough reputation to "comment" on this post. This webpage might be of interest to you. The development team describses here, although quite shortly, the main implications and solutions to warnings in STAN. Reaching the Maximum treedetph is, according to their explanation, far less an issue than divergent ...


5

The posterior predictive distribution of a logit model involves four conceptual steps: Draw from the posterior distribution of the parameters Create linear predictors using those parameter draws Transform those linear predictors into probabilities via the inverse link function Draw from a Bernoulli (or Binomial with size = 1) distribution with those ...


5

In BUGS/JAGS, the order in which statements are written does not matter. In Stan statements execute in the order in which they are written (see Stan 2.2.0 Reference Manual, pg. 405). Thus your last statement is in the wrong place: lambda is sampled from a gamma distribution, but that happens after the previous statements, so it's sampled 'in the air'. ...


5

First stan_lm only accepts one type of prior: Must be a call to R2 with its location argument specified or NULL, which would indicate a standard uniform prior for the R^2. If you want to specify any of the other types of priors, you'll need to switch to stan_glm. With the default family (gaussian), it will be equivalent to stan_lm. Second, you've tried to ...


5

$y \mid p,\lambda$ is Poisson! Your marginalization, or at least the end result, is correct. The form you have obtained for the distribution is the probability mass function of a Poisson distribution -- just write $p^y\lambda^y$ as $(\lambda\,p)^y$ and behold. That is, \begin{equation} y \mid p, \lambda \sim \mathrm{Poisson}(\lambda\,p). \end{equation} ...


5

Most software such as Stan, WingBUGS, SAS etc. requires you to provide an analytic form for the prior instead of MCMC samples. Possible ways around it are to refit the model with all data or to approximate the posterior with e.g. some mixture distribution (e.g. of bivariate normals to $\mu $ and $\log \sigma$ - e.g. in R using e.g. the mclust package or in ...


5

I don't know about pyStan specifically, but it's straight forward to compute an HDI from an MCMC sample, if you assume that the underlying distribution is unimodal. Basically, to compute the 95% HDI, you sort sample, then check all the 95% intervals, and select the shortest. Functions for computing HDI's are explained in Chapter 25 of Doing Bayesian Data ...


5

Each value of y_miss can either be 0 or 1, so you need to marginalize over them with a statement such as for (i in 1:n_miss) { real eta = a + beta * x_miss[i]; target += log_mix(inv_logit(eta), bernoulli_logit_lpmf(1 | eta), bernoulli_logit_lpmf(0 | eta)); } where x_miss is a vector of observations on the predictor ...


Only top voted, non community-wiki answers of a minimum length are eligible