# Estimating expectations of marginalised discrete variable

**URL:** <https://discourse.mc-stan.org/t/estimating-expectations-of-marginalised-discrete-variable/36934>\
**Category:** Modeling\
**Created:** [October 15, 2024, 2:57pm UTC](https://discourse.mc-stan.org/t/estimating-expectations-of-marginalised-discrete-variable/36934 "2024-10-15T14:57:34Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![jandraor](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/jandraor/32/7919_2.png) [@jandraor](https://discourse.mc-stan.org/u/jandraor)\
**Post date:** [October 15, 2024, 2:57pm UTC](https://discourse.mc-stan.org/t/estimating-expectations-of-marginalised-discrete-variable/36934/1 "2024-10-15T14:57:34Z")

</div>

Hi,

I’m following Stan’s manual ([Latent Discrete Parameters](https://mc-stan.org/docs/stan-users-guide/latent-discrete.html)) to estimate a discrete parameter. All chains converge and the method seems to recover the ground truth. At least visually, but I don’t know how to summarise it with a probability statement. Let’s assume that the parameter of interest is \tau, which can take values between 360 and 4920. I have evaluated the likelihood at discrete points of this interval in increments of 120.

This is the Stan code that generates the probabilities and draws. According to Stan’s manual, `sim_cp` should only be used for visual inspection, but for precise calculations (expectations), I should use `prob_cp`. I would like a statement that says “There is a 95% chance that the change point (cp) is between a & b”. Any pointers as to how to do that will be highly appreciated. I got lost with the notation in the section “estimating expectations”.

```no-highlight
generated quantities {
  array[n_i_y] int sim_cp;
  array[n_i_y] vector[n_cp] prob_cp;

  for(i in 1:n_i_y) {
    sim_cp[i] = cp[categorical_logit_rng(lp[i,])];
    prob_cp[i] = softmax(lp[i]);
  }
}

```

In R, for a given individual i, I have a data frame with three columns ([sample.csv](https://discourse.mc-stan.org/uploads/short-url/iTFJFXK5yp0FOBw0lVjZMlbhaYh.csv) (1.3 KB)), where `j` is the index of the change point, `prob_cp the` probability, and `.iter` the iteration.

Thanks.

---

<div class="post-metadata">

**Author:** ![Bob\_Carpenter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bob_carpenter/32/9230_2.png) [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)\
**Post date:** [October 21, 2024, 7:21pm UTC](https://discourse.mc-stan.org/t/estimating-expectations-of-marginalised-discrete-variable/36934/2 "2024-10-21T19:21:11Z")

</div>

> [@jandraor](#):
>
> According to Stan’s manual, `sim_cp` should only be used for visual inspection

You also want to add randomness here if you are doing posterior predictive inference.

In general, suppose you have an unknown `theta` in Stan. You can define an indicator as follows that assigns a value 1 if \theta \in (a, b) and value 0 otherwise (Stan’s comparison functions return either 0 for false or 1 for true).

```stan
generated quantities {
  int<lower=0, upper=1> Pr_btw = a < theta && theta < b;
}

```

The posterior mean of this quantity is the probability:

\Pr[\theta \in (a, b) \mid y] = \int\_\Theta \textrm{I}\_{(a, b)}\!(\theta) \cdot p(\theta | y) \, \textrm{d} \theta \approx \frac{1}{M} \sum\_{m=1}^M \textrm{I}\_{(a, b)}(\theta^{(m)}),

where \theta^{(m)} \sim p(\theta \mid y) are the posterior draws and \textrm{I}\_A(u) = 1 if u \in A and 0 otherwise (i.e., it’s an indicator function for the interval).