22.1 Softmax regression

The usual generalizations of logistic regression are of two types: conditional logistic regression, which will be explored later in the chapter, and so-called multinomial logistic regression, which we explore in this section. I am not fond of the traditional name, however, because the model does not use the logistic function per se, so the name “logistic” is a misnomer, and all the models in this chapter describe multinomial data, so the name “multinomial” is not informative. Instead, the key descriptor of the model is its inverse-link function, which is the softmax function (which will be defined below). Therefore, I refer to the method as softmax regression instead of multinomial logistic regression.

In the previous chapter, we used the logistic function to go from a linear combination of predictors to the probability of outcome 1 relative to outcome 0. In this chapter, we want to generalize to multiple categorical outcomes. The generalization of the logistic function requires a bit of mathematical development, but it's really just repeated algebraic manipulation of exponential functions, so don't be deterred. To keep the notation reasonably simple, we will suppose that we have a single metric predictor, x. All the derivations below easily generalize to multiple predictors. The underlying linear propensity of outcome k is denoted

${\text{λ}}_k={\beta }_{0,k}+{\beta }_{1,k}x$

The subscripts k indicate that there is an equation like Equation 22.1 for every outcome category. We call the set of possible outcomes S. Now a novelty: The probability of outcome k is given by the softmax function:

${\varphi }_k={\text{softmax}}_S\left(\left\{{\text{λ}}_k\right\}\right)=\frac{\exp \text{?}\left({\text{λ}}_k\right)}{{\Sigma }_{c\in S}\text{?}\exp \text{?}\left({\text{λ}}_c\right)}$

In words, Equation 22.2 says that the probability of outcome k is the exponentiated linear propensity of outcome k relative to the sum of exponentiated linear propensities across all outcomes in the set S. You may be wondering, Why exponentiate? Intuitively, we have to go from propensities that can have negative values to probabilities that can only have non-negative values, and we have to preserve order. The exponential function satisfies that need.

The softmax function is used in many applications as a way of mapping several real-valued variables to order-preserving outcome probabilities (e.g., Bishop, 2006). It is called the softmax function because when it is given another parameter, called the gain γ, that amplifies all the inputs, then the softmax assigns the maximum input with nearly 100% probability when the gain is large:

$\text{As}\text{?}\text{?}\text{?}\gamma \to \infty \text{,}\text{?}\text{?}\frac{\exp \left(\gamma {\text{λ}}_k\right)}{{\Sigma }_{c\in S}\text{?}\exp \text{?}\left(\gamma {\text{λ}}_c\right)}\to \{\begin{array}{l}1\text{?}\text{if}\text{?}{\text{λ}}_k=\text{max}\left(\left\{{\text{λ}}_c\right\}\right)\hfill \\ 0\text{}\text{?}\text{?}\text{}\text{otherwise}\hfill \end{array}\text{?}$

The softmax formulation can be useful for applications that use the derivative (i.e., gradient) to find optimal input values, because gradient ascent requires smooth differentiable functions (cf., Kruschke & Movellan, 1991).

It turns out that there are indeterminacies in the system of Equations 22.1 and 22.2. We can add a constant C₀ to every β_0,k, and add a constant C₁ to every β_{1, k}, and get exactly the same probabilities of responding with each category:

$\begin{array}{l}\frac{\exp \text{?}\left(({\beta }_{0,k}+C_0)+({\beta }_{1,k}+C_1)x\right)}{{\Sigma }_{c\in S}\exp \text{?}\left(({\beta }_{0,c}+C_0)+({\beta }_{1,c}+C_1)x\right)}\hfill \\ \text{?}\text{?}\text{?}=\frac{\exp \text{?}\left((C_0+C_{1}x)+{\beta }_{0,k}+{\beta }_{1,k}x\right)}{{\Sigma }_{c\in S}\text{?}\exp \text{?}\left((C_0+C_{1}x)+{\beta }_{0,c}+{\beta }_{1,c}x\right)}\hfill \\ \text{?}\text{?}\text{?}=\frac{\exp \text{?}(C_0+C_{1}x)\text{?}\exp \text{?}({\beta }_{0,k}+{\beta }_{1,k}x)}{{\Sigma }_{c\in S}\text{?}\exp \text{?}(C_0+C_{1}x)\text{?}\exp \text{?}({\beta }_{0,c}+{\beta }_{1,c}x)}\hfill \\ \text{?}\text{?}\text{?}=\frac{\exp \left({\beta }_{0,k}+{\beta }_{1,k}x\right)}{{\displaystyle {\sum }_{c\in S}\exp \left({\beta }_{0,c}+{\beta }_{1,c}x\right)}}\hfill \\ \text{?}\text{?}\text{?}={\varphi }_k\hfill \end{array}$

Therefore, we can set the baseline and slope for one of the response categories to arbitrary convenient constants. We will set the constants of one response category, called the reference category r, to zero: β_{0, r} =0 and β_{1, r} =0.

Because of the indeterminacy in the regression coefficients, we can interpret the regression coefficients only relative to the reference category. $\mathrm{R}$ecall from Section 21.2.1 (p. 629) that the regression coefficients in logistic regression can be conceived in terms of log odds of outcome 1 relative to outcome 0. In the present application, the regression coefficients can be conceived in terms of the log odds of each outcome relative to the reference outcome:

$\begin{array}{lll}\log \text{?}\left(\frac{{\varphi }_k}{{\varphi }_r}\right)\hfill & =\hfill & \log \text{?}\left(\frac{\exp \text{?}({\beta }_{0,k}+{\beta }_{1,k}x)}{\exp \text{?}({\beta }_{0,r}+{\beta }_{1,r}x)}\right)\hfill \\ \hfill & =\hfill & \log \text{?}\left(\frac{\exp \text{?}({\beta }_{0,k}+{\beta }_{1,k}x)}{\exp \text{?}(0+0x)}\right)\hfill \\ \hfill & =\hfill & {\beta }_{0,k}+{\beta }_{1,k}x\hfill \end{array}$

In other words, the regression coefficient β_{1, k} is the increase in log odds of outcome k relative to the reference outcome r for a one-unit increase in x.

Figure 22.1 shows two examples of data generated from the multinomial logistic model. The examples use two predictors (x₁ and x₂) to better illustrate the variety of patterns that can be produced by the model. The regression coefficients are displayed in the titles of the panels within the figure. The regression coefficients are large on the scale of x so that the transitions across outcomes are visually crisp. In realistic data, the regression coefficients are usually of much smaller magnitude and the outcomes of all types extensively overlap.

For the examples in Figure 22.1, the reference outcome was chosen to be 1, so all the regression coefficients describe the log odds of the other outcomes relative to outcome 1. Consider outcome 2, which has log odds specified by the regression coefficients for λ₂. If, for λ₂, the regression coefficient on x₁ is positive, then the probability of outcome 2, relative to outcome 1, goes up as x₁ increases. And, still for λ₂, if the regression coefficient on x₂ is positive, then the probability of outcome 2, relative to outcome 1, goes up as x₂ increases. The baseline, β₀ determines how large x₁ and x₂ need to become for the probability of outcome 2 to exceed the probability of outcome 1: The larger the baseline, the smaller x₁ and x₂ need to be.

You can see these trends occur in the two panels of Figure 22.1. In each panel, locate the region in which the reference outcome (i.e., 1) occurs. If λ₂'s coefficient on x₁ is positive, then the region in which outcome 2 occurs will tend to be on the right of the region where the reference outcome occurs. If λ₂'s coefficient on x₁ is negative, then the region in which outcome 2 occurs will tend to be on the left of the region where the reference outcome occurs. Analogously for x₂: If λ₂'s coefficient on x₂ is positive, then the region in which outcome 2 occurs will tend to be above the region where the reference outcome occurs. If λ₂'s coefficient on x₂ is negative, then the region in which outcome 2 occurs will tend to be below the region where the reference outcome occurs. These trends apply to the other outcomes as well.

Another general conclusion to take away from Figure 22.1 is that the regions of the different outcomes do not necessarily have piecewise linear boundaries. Later in the chapter we will consider a different sort of model that always produces outcome regions that have piecewise linear boundaries.

22.2 Conditional logistic regression

Softmax regression conceives of each outcome as an independent change in log odds from the reference outcome, and a special case of that is dichotomous logistic regression. But we can generalize logistic regression another way, which may better capture some patterns of data. The idea of this generalization is that we divide the set of outcomes into a hierarchy of two-set divisions, and use a logistic to describe the probability of each branch of the two-set divisions. The underlying propensity to respond with any outcome in the subset of outcomes S* relative to (i.e., conditional on) a larger set S is denoted

${\text{λ}}_{S*|S}={\beta }_{0,S*|S}+{\beta }_{1,S*|S}\text{?}x$

and the conditional response probability is

${\varphi }_{S*|S}=\text{logistic}\left({\text{λ}}_{S*|S}\right)$

The thrust of Equations 22.7 and 22.8 is that the regression coefficients refer to the conditional probability of outcomes for the designated subsets, not necessarily to a single outcome among the full set of outcomes.

Figure 22.2 shows two examples of hierarchies of binary divisions of four outcomes. In each example, the full set of (four) outcomes appears at the top of the hierarchy, and the binary divisions progress downward. Each branch is labeled with its conditional probability. For example, the top left branch indicates the occurrence of outcome 1 given the full set, and the conditional probability of this occurrence is labeled ϕ_{{1}|{1,2,3,4}}. This conditional probability will be modeled by a logistic function of the predictors. We will now explore detailed numerical examples of the two hierarchies.

In the left hierarchy of Figure 22.2, we opt to split the outcomes into a hierarchy of binary divisions as follows:

• 1 versus 2, 3, or 4

• 2 versus 3 or 4

• 3 versus 4

At each level in the hierarchy, the conditional probability of the options is described by a logistic function of a linear combination of the predictors. For our example, we assume there are two metric predictors. The linear combinations of predictors are denoted

$\begin{array}{ccc}\hfill {\lambda }_{\left\{1\right\}\left|\right\{1,2,3,4\}}& =& {\beta }_{0,\left\{1\right\}\left|\right\{1,2,3,4\}}+{\beta }_{1,\left\{1\right\}\left|\right\{1,2,3,4\}}{x}_1+{\beta }_{2,\left\{1\right\}\left|\right\{1,2,3,4\}}{x}_2\hfill \\ \hfill {\lambda }_{\left\{2\right\}\left|\right\{2,3,4\}}& =& {\beta }_{0,\left\{2\right\}\left|\right\{2,3,4\}}+{\beta }_{1,\left\{2\right\}\left|\right\{2,3,4\}}{x}_1+{\beta }_{2,\left\{2\right\}\left|\right\{2,3,4\}}{x}_2\hfill \\ \hfill {\lambda }_{\left\{3\right\}\left|\right\{3,4\}}& =& {\beta }_{0,\left\{3\right\}\left|\right\{3,4\}}+{\beta }_{1,\left\{3\right\}\left|\right\{3,4\}}{x}_1+{\beta }_{2,\left\{3\right\}\left|\right\{3,4\}}{x}_2\hfill \end{array}$

and the conditional probabilities of the outcome sets are simply the logistic function applied to each of the λ values:

$\begin{array}{ccc}\hfill {\varphi }_{\left\{1\right\}\left|\right\{1,2,3,4\}}& =& \text{logistic}({\lambda }_{\left\{1\right\}\left|\right\{1,2,3,4\}})\hfill \\ \hfill {\varphi }_{\left\{2\right\}\left|\right\{2,3,4\}}& =& \text{logistic}({\lambda }_{\left\{2\right\}\left|\right\{2,3,4\}})\hfill \\ \hfill {\varphi }_{\left\{3\right\}\left|\right\{3,4\}}& =& \text{logistic}({\lambda }_{\left\{3\right\}\left|\right\{3,4\}})\hfill \end{array}$

The ϕ values, above, should be thought of as conditional probabilities, as marked on the arrows in Figure 22.2. For example, ϕ_{2}|{2,3,4} is the probability of outcome 2 given that outcomes 2, 3, or 4 occurred. And the complement of that probability, 1 – ϕ_{2}|{2,3,4} is the probability of outcomes 3 or 4 given that outcomes 2, 3, or 4 occurred.

Finally, the most important equations that actually determine the relations between the expressions above, are the equations that specify the probabilities of the individual outcomes as the appropriate combinations of the conditional probabilities. These probabilities of the individual outcomes are determined simply by multiplying the conditional probabilities on the branches of the hierarchy in Figure 22.2 that lead to the individual outcomes. Thus,

$\begin{array}{l}{\varphi }_1={\varphi }_{\left\{1\right\}|\left\{1,2,3,4\right\}}\\ {\varphi }_2={\varphi }_{\left\{2\right\}|\left\{2,3,4\right\}}\cdot \left(1-{\varphi }_{\left\{1\right\}|\left\{1,2,3,4\right\}}\right)\\ {\varphi }_3={\varphi }_{\left\{3\right\}|\left\{3,4\right\}}\cdot \left(1-{\varphi }_{\left\{2\right\}|\left\{2,3,4\right\}}\right)\cdot \left(1-{\varphi }_{\left\{1\right\}|\left\{1,2,3,4\right\}}\right)\\ {\varphi }_4=\left(1-{\varphi }_{\left\{3\right\}|\left\{3,4\right\}}\right)\cdot \left(1-{\varphi }_{\left\{2\right\}|\left\{2,3,4\right\}}\right)\cdot \left(1-{\varphi }_{\left\{1\right\}|\left\{1,2,3,4\right\}}\right)\end{array}$

Notice that the sum of the outcome probabilities is indeed 1 as it should be; that is, ϕ₁ + ϕ₂ + ϕ₃ + ϕ₄ = 1. This result is easiest to verify in this case by summing in the order ((ϕ₄+ϕ₃)+ϕ₂)+ϕ₁.

The left panel of Figure 22.3 shows an example of data generated from the conditional logistic model of Equations 22.9–22.11. The x₁ and x₂ values were randomly generated from uniform distributions. The regression coefficients are displayed in the title of the panel within the figure. The regression coefficients are large compared to the scale of x so that the transitions across outcomes are visually crisp. In realistic data, the regression coefficients are usually of much smaller magnitude and outcomes of all types extensively overlap. You can see in the left panel of Figure 22.3 that there is a linear separation between the region of 1 outcomes and the region of 2, 3, or 4 outcomes. Then, within the region of 2, 3, or 4 outcomes, there is a linear separation between the region of 2 outcomes and the region of 3 or 4 outcomes. Finally, within the region of 3 or 4 outcomes, there is a linear separation between the 3's and 4's. Of course, there is probabilistic noise around the linear separations; the underlying linear separations are plotted where the logistic probability is 50%.

We now consider another example of conditional logistic regression, again involving four outcomes and two predictors, but with a different parsing of the outcomes. The hierarchy of binary partitions of the outcomes is illustrated in the right side of Figure 22.2. In this case, we opt to split the outcomes into the following hierarchy of binary choices:

• 1 or 2 versus 3 or 4

• 1 versus 2

• 3 versus 4

At each level in the hierarchy, the conditional probability of the options is described by a logistic function of a linear combination of the predictors. For the choices above, the linear combinations of predictors are denoted

$\begin{array}{ccc}\hfill {\lambda }_{\{1,2\}\left|\right\{1,2,3,4\}}& =& {\beta }_{0,\{1,2\}\left|\right\{1,2,3,4\}}+{\beta }_{1,\{1,2\}\left|\right\{1,2,3,4\}}{x}_1+{\beta }_{2,\{1,2\}\left|\right\{1,2,3,4\}}{x}_2\hfill \\ \hfill {\lambda }_{\left\{1\right\}\left|\right\{1,2\}}& =& {\beta }_{0,\left\{1\right\}\left|\right\{1,2\}}+{\beta }_{1,\left\{1\right\}\left|\right\{1,2\}}{x}_1+{\beta }_{2,\left\{1\right\}\left|\right\{1,2\}}{x}_2\hfill \\ \hfill {\lambda }_{\left\{3\right\}\left|\right\{3,4\}}& =& {\beta }_{0,\left\{3\right\}\left|\right\{3,4\}}+{\beta }_{1,\left\{3\right\}\left|\right\{3,4\}}{x}_1+{\beta }_{2,\left\{3\right\}\left|\right\{3,4\}}{x}_2\hfill \end{array}$

and the conditional probabilities of the outcome sets are simply the logistic function applied to each of the λ values:

$\begin{array}{ccc}\hfill {\varphi }_{\{1,2\}\left|\right\{1,2,3,4\}}& =& \text{logistic}({\lambda }_{\{1,2\}\left|\right\{1,2,3,4\}})\hfill \\ \hfill {\varphi }_{\left\{1\right\}\left|\right\{1,2\}}& =& \text{logistic}({\lambda }_{\left\{1\right\}\left|\right\{1,2\}})\hfill \\ \hfill {\varphi }_{\left\{3\right\}\left|\right\{3,4\}}& =& \text{logistic}({\lambda }_{\left\{3\right\}\left|\right\{3,4\}})\hfill \end{array}$

Finally, the most important equations, that actually determine the relations between the expressions above, are the equations that specify the probabilities of the individual outcomes as the appropriate combinations of the conditional probabilities, which can be gleaned from the conditional probabilities on the branches of the hierarchy in the right side of Figure 22.2. Thus:

$\begin{array}{l}{\varphi }_1={\varphi }_{\left\{1\right\}|\left\{1,2\right\}}\cdot {\varphi }_{\left\{1,2\right\}|\left\{1,2,3,4\right\}}\\ {\varphi }_2=\left(1-{\varphi }_{\left\{1\right\}|\left\{1,2\right\}}\right)\cdot {\varphi }_{\left\{1,2\right\}|\left\{1,2,3,4\right\}}\\ {\varphi }_3={\varphi }_{\left\{3\right\}|\left\{3,4\right\}}\cdot \left(1-{\varphi }_{\left\{1,2\right\}|\left\{1,2,3,4\right\}}\right)\\ {\varphi }_4=\left(1-{\varphi }_{\left\{3\right\}|\left\{3,4\right\}}\right)\cdot \left(1-{\varphi }_{\left\{1,2\right\}|\left\{1,2,3,4\right\}}\right)\end{array}$

Notice that the sum of the outcome probabilities is indeed 1 as it should be; that is, ϕ₁ + ϕ₂ + ϕ₃ + ϕ₄ = 1. This result is easiest to verify in this case by summing in the order (ϕ₁ + ϕ₂) + (ϕ₃ + ϕ₄). The only structural difference between this example and the previous example is the difference between the forms of Equations 22.11 and 22.14. Aside from that, both examples involve coefficients on predictors where the meaning of the coefficients is determined by how their values are combined in Equations 22.11 and 22.14. The coefficients were given mnemonic subscripts in anticipation of the particular combinations Equations 22.11 and 22.14.

The right panel of Figure 22.3 shows an example of data generated from the conditional logistic model of Equations 22.11–22.14. The regression coefficients are displayed in the title of the panel within the figure. You can see in the right panel of Figure 22.3 that there is a linear separation between the region of 1 or 2 outcomes and the region of 3 or 4 outcomes. Then, within the region of 1 or 2 outcomes, there is a linear separation between the region of 1 outcomes and the region of 2 outcomes. Finally, within the region of 3 or 4 outcomes, there is a linear separation between the 3's and 4's. Of course, there is probabilistic noise around the linear separations; the linear separations are plotted as dashed lines where the logistic probability is 50%.

In general, conditional logistic regression requires that there is a linear division between two subsets of the outcomes, and then within each of those subsets there is a linear division of smaller subsets, and so on. This sort of linear division is not required of the softmax regression model, examples of which we saw in Figure 22.1. You can see that the outcome regions of the softmax regression in Figure 22.1 do not appear to have the hierarchical linear divisions required of conditional logistic regression. $\mathrm{R}$eal data can be extremely noisy, and there can be multiple predictors, so it can be challenging or impossible to visually ascertain which sort of model is most appropriate. The choice of model is driven primarily by theoretical meaningfulness.

22.4 Generalizations and variations of the models

The goal of this chapter is to introduce the concepts and methods of softmax and conditional logistic regression, not to provide an exhaustive suite of programs for all applications. Fortunately, it is usually straight forward in principle to program in JAGS or Stan whatever model you may need. In particular, from the examples given in this chapter, it should be easy to implement any softmax or conditional logistic model on metric predictors.

Extreme outliers can affect the parameter estimates of softmax and conditional logistic regression. Fortunately it is easy to generalize and implement the robust modeling approach that was used for dichotomous logistic regression in Section 21.3 (p. 635). The predicted probabilities of the softmax or conditional logistic model are mixed with a “guessing” probability as in Equation 21.2 (p. 635), with the guessing probability being 1 over the number of outcomes.

Variable selection can be easily implemented. Just as predictors in linear regression or logistic regression can be given inclusions parameters, so can predictors in softmax or conditional logistic regression. The method is implemented just as was demonstrated in Section 18.4 (p. 536), and the same caveats and cautions still apply, as were explained throughout that section including subsection 18.4.1 regarding the influence of the priors on the regression coefficients.

The model can have nominal predictors instead of or in addition to metric predictors. For inspiration, consult the model diagram in Figure 21.12 (p. 642). A thorough development of this application would involve discussion of multinomial and Dirichlet distributions, which are generalizations of binomial and beta distributions. But that would take me beyond the intended scope of this chapter.

22.5 Exercises

Look for more exercises at https://sites.google.com/site/doingbayesiandataanalysis/