Question on cutpoint prior for stan_polr

[Kevin Van Horn]

I'm reading through the vignette "Estimating Ordinal Regression Models with rstanarm," and in the discussion of the cutpoint prior I'm getting confused as to what we are conditioning on when using the distribution for y^*.

It says that we specify a Dirichlet prior on \pi, where

  \pi_j = Pr(y = j | \bar{x})

I understand this to mean the conditional probability that y = j when x is set to its mean in the data set.

Later on it says that the j-th cutpoint is given by F^{-1}(p), where

  p = SUM(i: 1 <= i <= j: \pi_j))

and F is the cumulative probability distribution for y^*. Which distribution?

1. The distribution p(y^* | R)?

2. The distribution p(y^* | R, x = \bar{x})?

The definition of \pi implies that (2) is the distribution used to define the cutpoints... but then the vignette states that the scale parameter is \sigma_{y^*}, which implies that (1) is the defining distribution.

[Ben Goodrich]

On Monday, February 15, 2016 at 9:55:40 PM UTC-5, Kevin Van Horn wrote:

I'm reading through the vignette "Estimating Ordinal Regression Models with rstanarm," and in the discussion of the cutpoint prior I'm getting confused as to what we are conditioning on when using the distribution for y^*.

It says that we specify a Dirichlet prior on \pi, where

  \pi_j = Pr(y = j | \bar{x})

I understand this to mean the conditional probability that y = j when x is set to its mean in the data set.

Correct.
 

Later on it says that the j-th cutpoint is given by F^{-1}(p), where

  p = SUM(i: 1 <= i <= j: \pi_j))

and F is the cumulative probability distribution for y^*. Which distribution?

This depends on whether you are estimating an ordered logit model, an ordered probit model, etc. In the case of an ordered logit model, F() is the scaled logistic CDF.

This is one of the few cases where the relevant code is actually pretty readable. It is the make_cutpoints() function at

https://github.com/stan-dev/rstanarm/blob/master/exec/polr.stan

1. The distribution p(y^* | R)?

2. The distribution p(y^* | R, x = \bar{x})?

The definition of \pi implies that (2) is the distribution used to define the cutpoints... but then the vignette states that the scale parameter is \sigma_{y^*}, which implies that (1) is the defining distribution.

You can choose the Dirichlet prior on the simplex vector whose typical element is \pi_j without conditioning on the R matrix. This prior also determines the distribution of the cumulative sum of these probabilities, and then you can use the scaled inverse CDF to figure out what set of cutpoints is consistent with those cumulative probabilities.

Ben
 

[Kevin Van Horn]

I took a look at polr.stan, and it seems inconsistent with the claim that

  \pi_j = Pr(y = j | x = \bar{x})

The above equation is true for the cut-points computed in make_cutpoints() only if

1. \bar{x} = 0, so that E[y^*] = 0, and

2. scale = 1.

[Ben Goodrich]

On Tuesday, February 16, 2016 at 11:13:43 PM UTC-5, Kevin Van Horn wrote:

I took a look at polr.stan, and it seems inconsistent with the claim that

  \pi_j = Pr(y = j | x = \bar{x})

The above equation is true for the cut-points computed in make_cutpoints() only if

1. \bar{x} = 0, so that E[y^*] = 0, and

2. scale = 1.

The columns of the X matrix are centered and there is no intercept, so the first condition holds. And there is a scale parameter that goes into the make_cutpoints() function that multiplies the unit-scale cutpoints, so the second condition holds.

Ben
 

[Kevin Van Horn]

Yes, I know that there is a scale parameter passed to make_cutpoints; that is what I was referring to when I said that the condition

  \pi_j = Pr(y = j | x = \bar{x})

holds only if scale = 1, that is, only if the scale parameter passed in to make_cutpoints() equals 1. The reason it only holds when scale = 1 is because the error scale itself is 1, as described in the vignette, and as an examination of the pw_polr() and CDF_polr() functions confirms.

But the scale parameter passed in to make_cutpoints() does NOT equal 1 -- the Stan model uses 1/sqrt(1 - R^2), a value that is equal to 1 only if R = 0, and is otherwise strictly greater than 1.

[Jonah Gabry]

On Wednesday, February 17, 2016 at 11:08:05 PM UTC-5, Kevin Van Horn wrote:

The reason it only holds when scale = 1 is because the error scale itself is 1

 

Hi Kevin, can you clarify what you mean by this? Indeed the error scale (sigma_epsilon) is fixed at 1 since it isn't identified. But I'm not sure why if the unit-scaled cut-points are subsequently scaled by something other than 1 it invalidates the claim about the definition of pi. Thanks, 

Jonah

 

[Kevin Van Horn]

Let's step through the math. Since x is centered, the condition x = \bar{x} become x = 0, and P(y = k | x = \bar{x}) becomes P(y = k | x = 0). Then

    P(y = k | x = 0)

    = P(cutpoint[k-1] <= y^* <= cutpoint[k] | x = 0)

    = P(cutpoint[k-1] <= epsilon <= cutpoint[k])    # epsilon is the error term
    = F(cutpoint[k]) - F(cutpoint[k-1])                    # F is the CDF for the error distribution

If we define cutpoint[j] for all j to be F^{-1}(sum(\pi[1:j])), then we can continue as follows:

    = F(F^{-1}(sum(\pi[1:k]))) - F(F^{-1}(sum(\pi[1:(k-1)]))
    = sum(\pi[1:k]) - sum(\pi[1:(k-1)])
    = \pi[k]

If, however, we define cutpoint[j] = scale * F^{-1}(sum(\pi[1:j])) for scale != 1, then the equality no longer holds: F(scale * F^{-1}(z)) != z in general when scale != 1.