1
$\begingroup$

If I simulate the success probability for a Bernoulli random variable $t_{binary}$ conditional on some regressor $x$:

N <- 100000
x <- rnorm(N, 0, 1)
b0 <- 0
p <- plogis(b0 + 0.3*x)

the average probability (mean(p)) will be 0.5 as long as the intercept b0 is zero. Changing the slope of $x$ to something other than $0.3$ has no impact on this average.

But if I simulate the probabilities of a 3-category $t_{multi}$:

p0 <- exp(0 + 0.1*x)
p1 <- exp(0 + 0.2*x)
p2 <- exp(0 + 0.3*x)
p <- cbind(p0,p1,p2)
p <- t(apply(p, MARGIN = 1, function(h) h / sum(h)))

the average probability is split evenly (.33,.33,.33) only when the slope of x is the same for the three categories or if I fiddle with the intercepts. Is this correct (is this a property of the multinomial logit) or am I just simulating it wrong? The formula I have in mind is:

$P(t_{multi}=c|x)=\frac{exp(\beta_{0c}+\beta_cx)}{\sum_{k=0}^2exp(\beta_{0k}+\beta_kx)}$

Thanks.

$\endgroup$

1 Answer 1

2
$\begingroup$

This is a consequence of Jensen's inequality. Because the exponential function is convex, the expectation $\mathbb{E} \{ \exp(\beta_0 + \beta_1 x) \}$ is not the same as $\exp\{ \mathbb{E} (\beta_0 + \beta_1 x) \}$. Instead, $\mathbb{E} \{ \exp(\beta_0 + \beta_1 x) \} = \exp(\beta_0 + \frac{\beta_1^2}{2} )$

For logistic regression, this doesn't matter because zeroes and ones are symmetric due to how the logistic function is defined. But in your multinomial case, you're using the softmax function to set the probabilities, and so using your example of 0.1, 0.2, and 0.3 for the coefficients:

$ p(X=1) = \frac{ \exp(0.1^2/2 ) }{\exp(0.1^2/2 ) + \exp(0.2^2/2 ) + \exp(0.3^2/2 )} = 0.327$

$\endgroup$
2
  • 1
    $\begingroup$ Thank you! I had to Google the rule you mention at the end of your first paragraph. In case anyone else is wondering: since $x$ is standard normal, $\beta_0+\beta_1x$ is a normal with mean $\beta_0$ and standard deviation $\beta_1$; this is all inside an exponential function, so we have to apply an $ln$ to get to the normal stuff, meaning that $exp(\beta_0+\beta_1x)$ has a log-normal distribution. Then $exp(\beta_0+\beta_1^2/2)$ is simply the expectation of that log-normal, or $exp(\mu+\sigma^2/2)$. $\endgroup$
    – suckrates
    Commented Sep 27, 2021 at 8:44
  • 2
    $\begingroup$ I'm glad to hear you looked up the expectation of the exponential! I have a slightly different perspective than what you found, but it points to the same result: $\mathbb{E} \exp(\beta x)$ is the definition of the moment generating function of $X$. So you can find how to evaluate it by looking up the moment generating function of a normal random variable. And the $\exp(\beta_0)$ part is just a constant so you leave it as-is. $\endgroup$
    – Wesley
    Commented Sep 27, 2021 at 21:37

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.