Stats for

Generalized and Multilevel Regression

Multinomial Regression

Multinomial regression predicts which of several unordered categories someone falls into.

Multinomial regression takes an outcome with three or more unordered categories and tests how the odds of each category, compared with a chosen baseline category, differ with one or more predictors.

Example

Do more anxious students cope with exam stress by avoiding it, and do women and men cope differently? We asked 400 students how they mostly cope with exam stress (solving the problem, seeking support or avoiding it) and measured their trait anxiety on a 20–80 scale.

Unordered Outcome, 3+ Categories + Continuous or Categorical Predictors → Multinomial Regression
Tests
Population change in the log-odds of each category against the baseline, per unit of each predictor
Observations
Independent
Parametric
Yes
Headline Statistic
\chi^2
Effect Size
Relative Risk Ratio (RRR)
Groups
1

The Basics

When to Use It

Use it when your outcome is one of three or more categories with no natural order, each person falls into exactly one, and you want to know how one or more predictors (continuous or categorical) relate to which category people end up in. Also called multinomial logistic (or logit) regression, or polytomous logistic regression. For example:

  • Do more anxious students cope with exam stress by avoiding it, and do women and men cope differently?
  • Does severity of depressive symptoms (PHQ-9) predict whether clients choose individual therapy, group therapy or medication only?
  • Do age and years of education predict whether older adults' main memory complaint is names, appointments or where they put things?

Vocabulary You Need

Definition Baseline Category

In a model for an outcome with several unordered categories, the one category the others are each compared with. The model fits one equation per non-baseline category, and every coefficient and odds ratio describes that category rather than the baseline. Which category is the baseline doesn't change the model's fit or its predicted probabilities, only which comparisons are printed. Choose one that makes the comparisons meaningful: a common category, a natural default or the one your hypotheses are about. Software differs in its default (R's nnet::multinom uses the first level, SPSS the last, Stata the most frequent), so always say which you used. Also called the reference category or base outcome. Glossary →

Definition Log-Odds (\text{logit}(p))

The natural log of the odds, ln(p ÷ (1 − p)), also called the logit. A probability is stuck between 0 and 1 and odds between 0 and infinity, but log-odds run from minus to plus infinity, with 0 at a probability of .50: negative below it, positive above it, and symmetric, so .20 and .80 are −1.39 and +1.39. That unlimited scale is why logistic models are built on it: each predictor adds a fixed amount to the log-odds. A coefficient is a change in log-odds; e raised to it is an odds ratio, and 1 ÷ (1 + e^(−log-odds)) turns log-odds back into a probability. In a multinomial model the log-odds compare a category with the baseline category, ln(P(j) ÷ P(base)), rather than with everything else. Glossary →

Definition Odds Ratio (\text{OR})

One group's odds of an outcome divided by another group's. 1 means the odds are the same, 2 means the first group's odds are twice as high, and 0.5 means half as high; 0.5 and 2 are equally strong, in opposite directions. It's the usual effect size for a 2 × 2 table of two separate groups and what logistic regression reports. When the same people answer twice, as in McNemar's test, the paired odds ratio is different: it compares the people who switched one way with those who switched the other. It compares odds, not chances: when an outcome is common, the odds ratio is further from 1 than the ratio of the two proportions. A single zero cell, with the other counts positive, makes it 0 or infinite; more than one zero can leave it undefined. Glossary →

Definition Likelihood-Ratio Test (\chi^2)

Compares two nested models, a smaller one and a larger one that adds some parameters: the drop in deviance, χ² = D_smaller − D_larger = 2(LL_larger − LL_smaller), is compared with a chi-square distribution whose degrees of freedom are the number of added parameters (Wilks, 1938). Against a model with no predictors it tests the whole model at once, like a regression's F test; dropping one predictor tests that predictor alone. It is also written G² or LR χ². The Wald z tests that software prints beside each coefficient answer similar questions from a single fit, and in large samples the two usually agree. For mixed models, fit both models by maximum likelihood, not REML, when they differ in their fixed effects. When the tested value is a variance of zero, which sits on the edge of what a variance can be, the chi-square p-value is too conservative. Glossary →

Assumptions

1. Independent Observations The people (or other units) in your study don't influence each other.
What It Means
Each person's data comes from a separate, unconnected source, so knowing one person's scores tells you nothing about anyone else's. In designs that measure people more than once, their own scores are meant to be related; it's different people who must be independent. Methods built for clustered data, such as multilevel models and GEE, also allow members of the same group to be related, as long as the analysis names the grouping; then it's the groups that must be independent.
Applies To
Observations
How to Check
By design, not by a test. Did each person contribute their own data, with no clusters such as roommates, couples, therapy groups or classrooms whose members influence each other? If there are clusters, does the analysis model every one of them?
How Much It Matters
A lot. If students come in clusters (study groups, classes, friends who answered together), their choices are more alike than the model expects, the standard errors come out too small and every p is too small.

If It's Violated

  • Alternative Test Multilevel Logistic Regression. A multilevel (mixed-effects) multinomial or logistic model, which estimates the same comparisons while allowing for the clustering.
  • Design Change Sample people independently, or keep one person per cluster.
2. Independence of Irrelevant Alternatives Your outcome categories are distinct choices, not near-copies of each other.
What It Means
The odds of choosing one category over another don't depend on which other categories exist. Remove one strategy from the menu and the people who would have chosen it should spread over the rest in proportion to their current shares. That's plausible when the categories are clearly different; it fails when two categories are close substitutes, since people leaving one would move mostly to its twin.
Applies To
Outcome
How to Check
By reasoning about the categories: are any two of them near-substitutes that people would swap between? Tests exist, such as Hausman and McFadden's, which refit the model without one category and compare the coefficients, but in simulations they often reject IIA for distinct categories and miss it for near-substitutes, so don't let one decide.
How Much It Matters
It depends on the categories. The model assumes the odds of any two categories don't depend on what other categories are on offer. With clearly distinct strategies (problem-solving, seeking support, avoiding) that's usually plausible. If two categories are near-substitutes, such as two kinds of avoidance, the model misstates how people move between them. Formal tests, such as Hausman and McFadden's (1984), are unreliable in practice: their error rates depend on the predictors' data structure, and they often reject IIA for clearly distinct categories and miss it for near-substitutes (Cheng & Long, 2007), so the argument has to come from the categories themselves.

If It's Violated

  • Design Change Define categories that are clearly different from each other, and merge near-duplicates before collecting or analyzing the data.
  • Proceed with Caveat State that the categories were chosen to be distinct and that the model assumes IIA; don't lean on an IIA test to prove it.
3. Linearity of the Logit Each continuous predictor relates to the log-odds of the outcome in a straight line.
What It Means
The model assumes that each one-unit rise in a continuous predictor changes the log-odds of the outcome by the same amount, whether it happens low or high on the predictor's range. On the probability scale that's an S-shaped curve, so the probabilities themselves aren't expected to rise in a straight line. A U shape (risk that is high at both ends), a threshold or a leveling off breaks it. Yes/no predictors meet it automatically, since they have only two values.
Applies To
Relationship
How to Check
Add a bend and see whether it helps: a squared term, or a restricted cubic spline, for each continuous predictor, compared with the straight-line model by a likelihood-ratio test (Hosmer et al., 2013). The Box-Tidwell check adds x × ln(x) instead, for predictors with only positive values. A plot helps more than any p: cut the predictor into bands (such as deciles), and plot the log-odds of the outcome in each band against the band's average; the points should fall near a line.
How Much It Matters
A lot, for what the anxiety coefficients mean. Each equation assumes the log-odds of that strategy against problem-solving change by the same amount for every extra anxiety point. If avoidance only climbs at very high anxiety, one coefficient averages over that bend. Check each equation separately; the bend can differ between them. In the exam-coping study, adding anxiety² (centered) to both equations changes nothing, χ²(2) = 0.01, p = .996.

If It's Violated

  • Correction Add a squared term (anxiety², after centering anxiety) or a spline, and compare the models with a likelihood-ratio test.
  • Transformation Group the predictor into a few bands only if theory calls for it; it throws away information and the cut points change the answer.
4. No Severe Multicollinearity Your predictors aren't so closely related that their effects can't be told apart.
What It Means
Each coefficient is estimated from the part of its predictor that doesn't overlap with the others. When predictors are highly correlated, that part is small, so the coefficients become unstable and their standard errors large, even though the model as a whole may predict well. Perfect overlap, such as including a total score and all of its parts, makes the model impossible to estimate.
Applies To
Predictors
How to Check
Look at the correlations among the predictors and at each predictor's variance inflation factor (VIF). A VIF near 1 means little overlap. Common cutoffs of 5 or 10 are rules of thumb, not tests (O'Brien, 2007); the real question is whether the coefficients you care about are still precise enough to interpret.
How Much It Matters
Only if predictors overlap strongly. It inflates the standard errors in every equation at once. Anxiety and gender correlate only .12 in the exam-coping study, so it's not a concern here.

If It's Violated

  • Design Change If two predictors measure the same thing, average them into one score or keep the one that fits your question.
  • Proceed with Caveat If you only care about the predictors jointly, report the likelihood-ratio test for the set and don't lean on the individual coefficients.

Limitations

This Is Not Your Test If…

What It Doesn't Tell You

  • Why. Nothing here shows that anxiety makes students avoid. Students weren't assigned anxiety levels, and avoiding past exams could itself raise anxiety, or something left out, like workload, could drive both.
  • Whether it matters. The model improves on knowing nothing by only about 4% (McFadden's pseudo-R² = .039). Predicted chances of avoiding are about 16% at anxiety 35 and 31–35% at 55, two standard deviations higher, but most of what separates students' strategies isn't in the model: picking each student's most likely strategy is right 48% of the time, against 39% by always guessing problem-solving.
  • Who else. The results describe students like these 400, over the anxiety range they reported, and only for these three categories. Redefining them, or adding one that's close to an existing category, can change the comparisons.
  • That gender doesn't matter for avoiding. Women's non-significant coefficient for avoiding (p = .306, at the same anxiety level) means the data can't tell; its interval runs from 0.78 to 2.24, anywhere from a modest drop in women's odds of avoiding to more than double them.

Common Mistakes

  • Don'tRead a relative risk ratio as how many times more likely a category is.

    InsteadSay it multiplies the odds of that category rather than the baseline (despite Stata's name, it's an odds ratio, not a ratio of risks), and show predicted probabilities for the actual chances.

  • Don'tRead a positive coefficient as that category becoming more likely as the predictor rises.

    InsteadRead it as that category gaining on the baseline. Its own chance can still fall if another category gains faster, so show predicted probabilities across the predictor's range.

  • Don'tInterpret a coefficient without naming the baseline category it is compared with.

    InsteadSay “rather than problem-solving” (or whatever your baseline is) every time, and report which baseline you chose and why.

  • Don'tJudge a predictor by one equation's z test.

    InsteadTest it across all the equations at once with a likelihood-ratio test, then look at the separate comparisons.

  • Don'tFit separate logistic regressions for each pair of categories.

    InsteadFit one multinomial model, so every comparison uses all the data and the probabilities add up to 1.

  • Don'tRely on an IIA test to show the assumption holds.

    InsteadArgue from what the categories mean; IIA tests often give contradictory answers.

  • Don'tRead an odds ratio of 3 as “three times as likely.”

    InsteadSay “three times the odds,” and give the two percentages so readers can see the actual chances.

  • Don'tDescribe a regression coefficient as the effect of raising the predictor.

    InsteadSay the prediction is that much higher per unit of the predictor, on the model's scale (the outcome itself, or its log-odds or log rate in a logistic, ordinal or Poisson model). Causal language needs a design that supports it, such as an experiment.

  • Don'tRead a pseudo-R² as the share of variance explained, or compare it with a linear regression's R².

    InsteadName which pseudo-R² you report, and compare it only with the same kind from similar models.

  • Don'tRead a bigger χ² as a stronger effect.

    InsteadReport an effect size suited to the design, such as Cramér's V, an odds ratio or a model's coefficients. χ² grows with sample size even when the effect doesn't.

  • Don'tRead a p-value as the chance that the null hypothesis is true.

    InsteadRead it as: if the null were true, how often data at least this extreme would turn up.

Result Decoder

Significance and Effect Size Guidelines

Benchmarks Significance (p Value)
  • p < .001

    Very strong evidence against the null hypothesis.

  • .001 ≤ p < .01

    Strong evidence against the null hypothesis.

  • .01 ≤ p < .05

    Statistically significant at α = .05. Moderate evidence against the null hypothesis.

  • .05 ≤ p < .10

    Not significant. Weak evidence; sometimes called a “trend.” Treat it as inconclusive.

  • p ≥ .10

    Not significant. Little evidence against the null, which is not the same as evidence of no effect.

These evidence labels are one common convention for describing p in graded language rather than as significant or not; the boundaries are rough guides, not rules, and other authors draw them differently. Muff et al. (2022) argue for this kind of language of evidence over a single cutoff. α is chosen before the analysis, and p never says how big or important an effect is.

More about statistical significance ↓

Benchmarks Effect Size (Relative Risk Ratio (RRR))
  • 0.60 < OR < 1.68 Below Small

    Below the odds ratio that matches d = 0.2 when 1% of the comparison group has the outcome.

  • 1.68 ≤ OR < 3.47 (0.29 < OR ≤ 0.60) Small

    About d = 0.2 when 1% of the comparison group has the outcome.

  • 3.47 ≤ OR < 6.71 (0.15 < OR ≤ 0.29) Medium

    About d = 0.5 when 1% of the comparison group has the outcome.

  • OR ≥ 6.71 (OR ≤ 0.15) Large

    About d = 0.8 or more when 1% of the comparison group has the outcome.

An RRR is an odds ratio for one category against the baseline, so Chen, Cohen and Chen's (2010) cutoffs for odds ratios are borrowed here. They were derived for a two-group comparison, so they fit a yes/no predictor such as gender. For a continuous predictor the RRR depends on the unit: judge it for a meaningful change (here 10 anxiety points, about 1 SD), not for one point. For an RRR below 1, use 1 ÷ RRR or the bands in brackets. The cutoffs also assume a rare outcome (1% of the comparison group). Here every category is common (24–39% of students), and a common outcome needs a smaller ratio for the same Cohen's d, so the labels understate effects. Treat them as rough at best, and judge size from the predicted probabilities (Table 4).

More about effect size ↓

How It's Reported

Equation Form

Example

Model: χ²(4, N = 400) = 33.41, p < .001, McFadden's pseudo-R2 = .039; anxiety, avoid vs problem-solve: b = 0.061, RRR = 1.063 per point, 95% CI [1.034, 1.093], z = 4.36, p < .001; anxiety across both comparisons: χ²(2) = 20.61, p < .001

Table 1. Trait Anxiety and Gender by Coping Strategy

Example
nWomen (%)Anxiety MAnxiety SD
Problem-Solve 15650.042.379.71
Seek Support 14770.145.489.68
Avoid 9759.848.339.78
All Students 40059.844.969.97

Note. N = 400. Trait anxiety is a 20–80 questionnaire total.

Table 2. Likelihood-Ratio Tests

Example
−2LL without\chi^2dfp
Model (anxiety and gender) 862.9433.414< .001
Trait anxiety 850.1420.612< .001
Gender 840.3110.782.005

Note. Each test compares the full model with one that leaves the term out: for the model, a model with intercepts only. −2LL is that smaller model's deviance; the full model's is 829.53. Each predictor's test covers both of its coefficients, one per comparison with problem-solving.

Table 3. Multinomial Logistic Regression Predicting Coping Strategy (Baseline: Problem-Solve)

Example
bSEzpRRR95% CI
Seek Support vs Problem-Solve
Intercept −1.7970.563−3.19.001
Trait anxiety 0.0290.0122.33.0201.029[1.005, 1.054]
Woman 0.7850.2443.22.0012.193[1.360, 3.536]
Avoid vs Problem-Solve
Intercept −3.4030.665−5.12< .001
Trait anxiety 0.0610.0144.36< .0011.063[1.034, 1.093]
Woman 0.2760.2701.02.3061.318[0.777, 2.237]

Note. b = coefficient in log-odds; RRR = relative risk ratio, eb, with its 95% CI. Anxiety is per questionnaire point; woman is women against men. McFadden's pseudo-R2 = .039, Nagelkerke's = .091.

Table 4. Predicted Probability of Each Strategy

Example
Problem-SolveSeek SupportAvoid
Women, anxiety 35 .42.42.16
Women, anxiety 55 .25.44.31
Men, anxiety 35 .58.26.16
Men, anxiety 55 .36.29.35

Note. From the model in Table 3, at trait anxiety of 35 and 55, about 1 SD below and above the mean of 44.96. Each row adds to 1.

Full Text (APA Style)

Example

Of 400 students, 39% coped with exam stress mainly by problem-solving, 37% by seeking support and 24% by avoiding it. Students who avoided reported more trait anxiety (M = 48.33, SD = 9.78) than those who problem-solved (M = 42.37, SD = 9.71). In a multinomial logistic regression with problem-solving as the baseline, among students of the same gender, each additional anxiety point went with 1.063 times the odds of avoiding rather than problem-solving, 95% CI [1.034, 1.093], z = 4.36, p < .001, or 1.84 times per 10 points; for men, the predicted chance of avoiding was 16% at an anxiety score of 35 and 35% at 55. Higher anxiety was also associated with higher odds of seeking support rather than problem-solving, though less strongly, RRR = 1.029, 95% CI [1.005, 1.054], z = 2.33, p = .020, and across both comparisons it was statistically significant, χ²(2) = 20.61, p < .001. At the same anxiety level, women's odds of seeking support rather than problem-solving were 2.19 times men's, 95% CI [1.36, 3.54], z = 3.22, p = .001, but their odds of avoiding were not reliably different, RRR = 1.32, 95% CI [0.78, 2.24], z = 1.02, p = .306; gender was significant across both comparisons, χ²(2) = 10.78, p = .005. Together, anxiety and gender predicted coping strategy, χ²(4, N = 400) = 33.41, p < .001, but modestly, McFadden's pseudo-R2 = .039.

  • Name the baseline category, and write every comparison against it: “the odds of avoiding rather than problem-solving.”
  • Give a coefficient table like Table 3, with b, SE, z, p and the RRR with its 95% CI for every comparison, and the likelihood-ratio tests like Table 2. Say what software calls the ratio: RRR (Stata), Exp(B) (SPSS), or exp(coef) in R.
  • Translate the key results into predicted probabilities at meaningful predictor values, like Table 4, so readers can see the actual chances.
  • Name the pseudo-R² you report. McFadden's, Cox and Snell's and Nagelkerke's differ a lot for the same model.

Reading Each Piece

The equation form above, piece by piece, in this order: how big, how certain, then how surprising.

RRR = 1.063 Students one point higher in trait anxiety have 1.063 times the odds of avoiding rather than problem-solving, among students of the same gender; 10 points higher, about one standard deviation, 1.84 times. The benchmarks above would call that small, but they assume a category far rarer than avoiding (24% of students), so they understate it. It's a ratio of odds, not of chances; the predicted probabilities show the chances. Benchmarks ↑
95% CI [1.034, 1.093] Per anxiety point, the ratio of the odds of avoiding rather than problem-solving is plausibly anywhere from 1.034 to 1.093, or 1.40 to 2.43 per 10 points. The whole range is above 1.
p̂avoid = .35 (man, anxiety 55) A man with an anxiety score of 55 has about a 35% predicted chance of mainly avoiding, against 16% for a man at 35. For women the same 20 points take it from 16% to 31% (Table 4).
McFadden's pseudo-R2 = .039 The model improves the log-likelihood by about 4% over knowing only how common each strategy is. That's modest: anxiety and gender predict the odds, but most of what separates students' strategies isn't in the model. McFadden's values run much lower than ordinary R².
χ²(4, N = 400) = 33.41 Knowing students' anxiety and gender cuts the deviance by 33.41, far more than the 9.49 that chance exceeds 5% of the time with 4 df. The 4 is the degrees of freedom, 2 predictors × 2 comparisons, and N = 400 is the number of students.
Anxiety: χ²(2) = 20.61, p < .001 Taking anxiety out of both equations makes the deviance worse by 20.61. If anxiety's coefficients were 0 in both comparisons, a rise at least that large would turn up less than once in 10,000 studies. Gender's 10.78 (p = .005) says it's related to strategy in at least one comparison; Table 3 shows the clear difference is in seeking support, and for avoiding the data can't tell.
z = 4.36, p < .001 Anxiety's coefficient for avoiding is 4.36 standard errors above zero. In the support comparison it's 2.33 (p = .020), smaller but still significant (Table 3).
p < .001 If anxiety and gender predicted nothing about coping, a drop in deviance at least this big would turn up less than once in 1,000 studies (about 1 in a million). It says the predictors help; it doesn't say how much. Benchmarks ↑

Calculations

The Running Example

Example

The exam-coping study: 400 students, each with one main coping strategy (problem-solve, seek support or avoid), a trait anxiety score (20–80) and gender. Problem-solving is the baseline.

Null HypothesisAlternative
Modelall four slope coefficients are 0 (the intercepts stay)at least one is not 0
Anxietyb_{\text{anx, support}} = b_{\text{anx, avoid}} = 0at least one is not 0
Genderb_{\text{woman, support}} = b_{\text{woman, avoid}} = 0at least one is not 0

Each coefficient compares one strategy with problem-solving, so a null of 0 means that predictor isn't related to the odds of that strategy rather than problem-solving. The model test asks whether anxiety and gender together predict strategy at all; each predictor's test asks whether it matters in either comparison. Single coefficients also get their own two-sided z tests.

What Your Data Should Look Like

One row per person: one column naming their category (text or a numeric code), and one column for each predictor. A categorical predictor such as gender becomes 0/1 indicators, one fewer than its levels. Check how many people fall in each category, overall and for each level of a categorical predictor; categories with very few people give unstable coefficients. If a category has no one (or everyone) at some level of a predictor, the model can't estimate that coefficient: watch for huge coefficients and standard errors, or a warning that the model didn't converge. Merge sparse categories if that makes sense, or use a penalized (Firth-type) model.

Example
participantcopingtrait_anxietygender
1seek_support38woman
2problem_solve44woman
3seek_support48woman
4avoid40man
5problem_solve30man
6seek_support44woman
… 394 more rows
Problem-SolveSeek SupportAvoidTotal
Women78 (33%)103 (43%)58 (24%)239
Men78 (48%)44 (27%)39 (24%)161
Total15614797400

Counts, with each row's percentages.

Simulated data.

Download all 400 rows (data.csv)

Napkin Math

The whole test, one idea at a time. Under each heading, a line says why the step matters, and the green note reads the equation in plain English with the exam-coping study's numbers. Collapse all to read just the headings: they tell the whole story on their own.

  1. 1. Count How Often Students Chose Each Strategy and Pick the One to Compare Against

    A multinomial model compares every category with one baseline category. Problem-solving is the most common strategy and the healthy default, so the other two are each asked about against it.

    \text{odds}_{j} = \dfrac{n_j}{n_{\text{base}}} =\dfrac{97}{156}

    \text{odds}_{\text{avoid}} = 0.62

    \text{odds}_{\text{support}} = 0.94

    \hat{p}_{\text{problem}} = .39

    156 students mostly problem-solved (.39), 147 sought support (.37) and 97 avoided (.24). n_j is the 97 who avoided and n_{\text{base}} the 156 who problem-solved, so the odds of avoiding rather than problem-solving are 97 ÷ 156 = 0.62; for seeking support they're 147 ÷ 156 = 0.94.

  2. 2. Write One Equation for Each Strategy against Problem-Solving

    Each non-baseline strategy gets its own equation predicting its log-odds against problem-solving from anxiety and gender. Software finds all six coefficients at once by maximum likelihood, the values that make the students' actual choices most probable.

    \ln\dfrac{P(Y = j)}{P(Y = \text{base})} = b_{0j} + b_{1j}x_1 + b_{2j}x_2 =-3.40 + 0.061\,\text{anxiety} + 0.28\,\text{woman}

    b_{0,\text{avoid}} = −3.40

    b_{1,\text{avoid}} = 0.061

    b_{2,\text{avoid}} = 0.28

    b_{0,\text{support}} = −1.80

    b_{1,\text{support}} = 0.029

    b_{2,\text{support}} = 0.79

    That's the avoid equation; the support one is -1.80 + 0.029 anxiety + 0.79 woman. There's no formula for the answer, so software starts from zeros and improves the guesses (Newton-Raphson) until they stop changing, here after 5 rounds. b_{1,\text{avoid}} = 0.061 says students one point higher in anxiety have log-odds of avoiding rather than problem-solving 0.061 higher, holding gender constant.

  3. 3. Turn Anxiety's Coefficient into How Much It Multiplies the Odds

    Log-odds are hard to picture. Raising e to the coefficient turns it into a ratio of odds, which software calls the relative risk ratio (RRR) or Exp(B).

    \text{RRR} = e^{\,b_{kj}} =e^{0.0612}

    = 1.063

    \text{RRR}_{10\text{ points}} = 1.84

    b_{kj}, predictor k's coefficient in category j's equation, is anxiety's 0.0612 in the avoid equation (step 2's b_{1,\text{avoid}}). e^{0.0612} = 1.063: students one point higher in anxiety have 1.063 times the odds of avoiding rather than problem-solving, holding gender constant. Over 10 points, about 1 SD, that's e^{10 \times 0.0612} = 1.84.

  4. 4. Give a Range of Plausible Values for That Ratio

    Build the interval on the log-odds scale, where the coefficient's estimates from repeated samples would be roughly normal, so 1.96 standard errors either side covers 95%. Then turn both ends back into ratios.

    \left[\, e^{\,b - 1.96\,SE_b},\; e^{\,b + 1.96\,SE_b} \,\right] =\left[e^{0.0612 - 1.96 \times 0.0140},\; e^{0.0612 + 1.96 \times 0.0140}\right]

    = [1.034, 1.093]

    SE_b = 0.0140

    SE_b is 0.0140, the standard error software reads off the same fitting. 1.96 × 0.0140 = 0.0275, so the log-scale range is 0.0337 to 0.0887, and e to those is 1.034 to 1.093 per point.

  5. 5. Translate the Model Back into Chances of Each Strategy

    The clearest way to show what the coefficients mean is the predicted probability of each strategy for a particular student. Plug their values into both equations, then share out the probability. η (eta) is each equation's value for this student: step 2's right-hand side with their scores plugged in. The baseline's η is 0, which is where the 1 comes from, and the sum runs over the non-baseline categories.

    P(Y = j) = \dfrac{e^{\eta_j}}{1 + \sum_{m \ne \text{base}} e^{\eta_m}}, \quad P(Y = \text{base}) = \dfrac{1}{1 + \sum_{m \ne \text{base}} e^{\eta_m}} =\dfrac{e^{-0.037}}{1 + e^{-0.215} + e^{-0.037}}

    = .35

    \hat{p}_{\text{avoid}}(35) = .16

    For a man with anxiety 55, \eta_{\text{avoid}} = −3.403 + 0.0612 × 55 = −0.037 and \eta_{\text{support}} = −1.797 + 0.02876 × 55 = −0.215. So P(avoid) = 0.964 ÷ (1 + 0.806 + 0.964) = .35. At anxiety 35 the same arithmetic gives .16.

  6. 6. Ask Whether Anxiety and Gender Together Beat Knowing Nothing

    Compare the model with one that knows only how common each strategy is. The drop in deviance is the likelihood-ratio test; its degrees of freedom are the four coefficients the smaller model sets to zero.

    \chi^2_{LR} = D_{\text{smaller}} - D_{\text{larger}} =862.94 - 829.53

    = 33.41

    D_{\text{null}} = 862.94

    D_{\text{model}} = 829.53

    p < .001

    D_{\text{smaller}} is 862.94, −2 times the log-likelihood of the intercepts-only model, and D_{\text{larger}} is 829.53 for the full model. 862.94 − 829.53 = 33.41 on 2 predictors × 2 comparisons = 4 df. The chance a χ² with 4 df reaches 33.41 is under .001.

  7. 7. Measure How Much Better the Model Predicts than Knowing Nothing

    An R²-like number for the whole model. McFadden's pseudo-R² compares the two log-likelihoods; it is not a share of variance explained.

    R^2_{\text{McF}} = 1 - \dfrac{LL_{\text{model}}}{LL_{\text{null}}} =1 - \dfrac{-414.76}{-431.47}

    = .039

    R^2_{\text{N}} = .091

    LL_{\text{model}} is −414.76 and LL_{\text{null}} −431.47. Their ratio is .961, so 1 − .961 = .039. SPSS also prints Nagelkerke's version, .091, for the same model.

  8. 8. Test Each Predictor across Both Comparisons at Once

    Anxiety has two coefficients, one per comparison with problem-solving. Dropping both and seeing how much the deviance rises tests anxiety as a whole, with 2 df.

    \chi^2_{LR} = D_{\text{smaller}} - D_{\text{larger}} =850.14 - 829.53

    = 20.61

    p_{\text{anxiety}} < .001

    \chi^2_{\text{gender}} = 10.78

    p_{\text{gender}} = .005

    Without anxiety, D_{\text{smaller}} = 850.14; with it, 829.53. 850.14 − 829.53 = 20.61 on 2 df, p < .001. Gender, tested the same way, gives 840.31 − 829.53 = 10.78, p = .005.

  9. 9. Check Whether Anxiety's Coefficient for Avoiding Could Be Zero

    Each single coefficient also gets its own test, the coefficient divided by its standard error, judged against the normal curve.

    z = \dfrac{b_{kj}}{SE_{b_{kj}}} =\dfrac{0.0612}{0.0140}

    = 4.36

    p < .001

    b_{kj} is anxiety's 0.0612 in the avoid equation and SE 0.0140, so z = 4.36. The chance a normal value lands at least 4.36 from zero, in either direction, is under .001.

Run It in Software

The code below reads the running example's data from data.csv; save it in the same folder you run the code from.

SPSS (version 30, as documented by IBM)
  1. Make a 0/1 indicator for gender first (Transform → Compute Variable: woman = (gender = 'woman')), so the coefficient compares women with men. As a Factor, SPSS would make the last category (woman) the reference instead.
  2. Analyze → Regression → Multinomial Logistic.
  3. Move coping into Dependent, click Reference Category and choose Custom: problem_solve. SPSS's default is the last category.
  4. Move trait_anxiety and woman into Covariate(s).
  5. Statistics: check Pseudo R-square, Model fitting information, Likelihood ratio tests and Parameter estimates (95% confidence interval). Click Continue, then OK.

Or paste this syntax:

COMPUTE woman = (gender = 'woman').
EXECUTE.
NOMREG coping (BASE='problem_solve' ORDER=ASCENDING) WITH trait_anxiety woman
  /CRITERIA CIN(95) DELTA(0) MXITER(100) MXSTEP(5) CHKSEP(20) LCONVERGE(0) PCONVERGE(0.000001) SINGULAR(0.00000001)
  /MODEL
  /STEPWISE=PIN(.05) POUT(0.1) MINEFFECT(0) RULE(SINGLE) ENTRYMETHOD(LR) REMOVALMETHOD(LR)
  /INTERCEPT=INCLUDE
  /PRINT=PARAMETER SUMMARY LRT CPS STEP MFI.

Reading the Output

Read the model first: Chi-Square and Sig. in Model Fitting Information and McFadden in Pseudo R-Square. Then each predictor's row in Likelihood Ratio Tests. Then, in Parameter Estimates, each comparison's Exp(B) (SPSS's name for the relative risk ratio) with its 95% Confidence Interval for Exp(B), then Wald and Sig. SPSS's Wald is z², judged against χ² with 1 df, so its p is the same as the z test's.

Model Fitting Information

Model −2 Log Likelihood Chi-Square df Sig.
Intercept Only 862.936
Final 829.527 33.409 4 <.001

Pseudo R-Square

Cox and Snell Nagelkerke McFadden
.080 .091 .039

Likelihood Ratio Tests

Effect −2 Log Likelihood of Reduced Model Chi-Square df Sig.
Intercept 860.015 30.488 2 <.001
trait_anxiety 850.137 20.610 2 <.001
woman 840.305 10.778 2 .005

Parameter Estimates (the reference category is problem_solve)

coping B Std. Error Wald df Sig. Exp(B) 95% CI Lower Upper
avoid Intercept -3.403 .665 26.213 1 <.001
trait_anxiety .061 .014 19.043 1 <.001 1.063 1.034 1.093
woman .276 .270 1.046 1 .306 1.318 .777 2.237
seek_support Intercept -1.797 .563 10.203 1 .001
trait_anxiety .029 .012 5.427 1 .020 1.029 1.005 1.054
woman .785 .244 10.371 1 .001 2.193 1.360 3.536

The Intercept rows are the log-odds for a man with anxiety 0, which no one scores, so they aren't interpreted. For the comparison between seeking support and avoiding, rerun with one of them as the reference category.

Output laid out from this page's computed results, not captured from a live run; labels and layout may differ slightly in your version.

R
library(nnet)
coping <- read.csv("data.csv")
# The first level is the baseline.
coping$coping <- factor(coping$coping, levels = c("problem_solve", "seek_support", "avoid"))
coping$gender <- factor(coping$gender, levels = c("man", "woman"))

fit <- multinom(coping ~ trait_anxiety + gender, data = coping, trace = FALSE)
summary(fit)

# z tests, relative risk ratios and their 95% CIs
z <- coef(fit) / summary(fit)$standard.errors
round(2 * pnorm(-abs(z)), 4)
round(exp(coef(fit)), 3)
round(exp(confint(fit)), 3)

# Likelihood-ratio tests: the whole model, then each predictor across both equations
null <- multinom(coping ~ 1, data = coping, trace = FALSE)
anova(null, fit)
car::Anova(fit)

# McFadden's pseudo R-squared
1 - as.numeric(logLik(fit)) / as.numeric(logLik(null))

# Predicted probabilities at anxiety 35 and 55 (about 1 SD below and above the mean)
new <- expand.grid(trait_anxiety = c(35, 55), gender = c("woman", "man"))
cbind(new, round(predict(fit, newdata = new, type = "probs"), 3))
Call:
multinom(formula = coping ~ trait_anxiety + gender, data = coping,
    trace = FALSE)

Coefficients:
             (Intercept) trait_anxiety genderwoman
seek_support   -1.797224    0.02876320   0.7851158
avoid          -3.403407    0.06121631   0.2760852

Std. Errors:
             (Intercept) trait_anxiety genderwoman
seek_support   0.5626432    0.01234733   0.2437944
avoid          0.6647523    0.01402820   0.2698874

Residual Deviance: 829.5268
AIC: 841.5268
             (Intercept) trait_anxiety genderwoman
seek_support      0.0014        0.0198      0.0013
avoid             0.0000        0.0000      0.3063
             (Intercept) trait_anxiety genderwoman
seek_support       0.166         1.029       2.193
avoid              0.033         1.063       1.318
, , seek_support

              2.5 % 97.5 %
(Intercept)   0.055  0.499
trait_anxiety 1.005  1.054
genderwoman   1.360  3.536

, , avoid

              2.5 % 97.5 %
(Intercept)   0.009  0.122
trait_anxiety 1.034  1.093
genderwoman   0.777  2.237

                   Model Resid. df Resid. Dev   Test    Df LR stat.
1                      1       798   862.9355           NA       NA
2 trait_anxiety + gender       794   829.5268 1 vs 2     4 33.40861
       Pr(Chi)
1           NA
2 9.851269e-07
Analysis of Deviance Table (Type II tests)

Response: coping
              LR Chisq Df Pr(>Chisq)
trait_anxiety   20.610  2  3.347e-05 ***
gender          10.778  2   0.004566 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
[1] 0.03871507
  trait_anxiety gender problem_solve seek_support avoid
1            35  woman         0.422        0.420 0.158
2            55  woman         0.248        0.438 0.315
3            35    man         0.576        0.261 0.163
4            55    man         0.361        0.291 0.348

multinom sets the first level's coefficients to zero, so set the factor levels to choose the baseline. It prints coefficients and standard errors but no tests; the lines after summary add the z tests, the RRRs (exp(coef)) and their intervals. nnet fits by an optimizer that stops a little early, so its coefficients differ from other programs' in the sixth or seventh digit. car::Anova gives each predictor's likelihood-ratio test across both equations.

Output from a real run: R 4.6.1 (nnet 7.3-20, car 3.1-5), run 2026-09-29.

Python
import numpy as np
import pandas as pd
import statsmodels.api as sm

coping = pd.read_csv("data.csv")
# The first category is the baseline.
strategy = pd.Categorical(coping["coping"], categories=["problem_solve", "seek_support", "avoid"])
X = pd.DataFrame({"trait_anxiety": coping["trait_anxiety"], "woman": (coping["gender"] == "woman").astype(int)})
X = sm.add_constant(X)

fit = sm.MNLogit(strategy.codes, X).fit()
print(fit.summary())

print(np.exp(fit.params))     # relative risk ratios
print(np.exp(fit.conf_int()))  # and their 95% CIs

# Likelihood-ratio test for anxiety across both equations
no_anxiety = sm.MNLogit(strategy.codes, X[["const", "woman"]]).fit(disp=0)
lr = 2 * (fit.llf - no_anxiety.llf)
print(lr)
Optimization terminated successfully.
         Current function value: 1.036909
         Iterations 5
                          MNLogit Regression Results
==============================================================================
Dep. Variable:                      y   No. Observations:                  400
Model:                        MNLogit   Df Residuals:                      394
Method:                           MLE   Df Model:                            4
Date:                Tue, 29 Sep 2026   Pseudo R-squ.:                 0.03872
Time:                        13:29:30   Log-Likelihood:                -414.76
converged:                       True   LL-Null:                       -431.47
Covariance Type:            nonrobust   LLR p-value:                 9.851e-07
=================================================================================
          y=1       coef    std err          z      P>|z|      [0.025      0.975]
---------------------------------------------------------------------------------
const            -1.7972      0.563     -3.194      0.001      -2.900      -0.694
trait_anxiety     0.0288      0.012      2.330      0.020       0.005       0.053
woman             0.7851      0.244      3.220      0.001       0.307       1.263
---------------------------------------------------------------------------------
          y=2       coef    std err          z      P>|z|      [0.025      0.975]
---------------------------------------------------------------------------------
const            -3.4034      0.665     -5.120      0.000      -4.706      -2.101
trait_anxiety     0.0612      0.014      4.364      0.000       0.034       0.089
woman             0.2761      0.270      1.023      0.306      -0.253       0.805
=================================================================================
                      0         1
const          0.165758  0.033259
trait_anxiety  1.029181  1.063129
woman          2.192660  1.317966
                    lower     upper
y
1 const          0.055024  0.499340
  trait_anxiety  1.004573  1.054391
  woman          1.359731  3.535814
2 const          0.009038  0.122392
  trait_anxiety  1.034297  1.092765
  woman          0.776561  2.236829
20.60972554785519

y=1 is seek support and y=2 is avoid, each against the baseline (code 0, problem-solve). The top block has the model's likelihood-ratio test (LLR p-value, from LL-Null and Log-Likelihood) and McFadden's pseudo-R² (Pseudo R-squ.). The interval in the table is for b; exponentiate it for the RRRs' intervals, as the last lines do.

Output from a real run: Python 3.9.6 (statsmodels 0.14.6, pandas 2.3.3), run 2026-09-29.

Impact

Effect Size

A multinomial model has effect sizes at three levels. Each coefficient's relative risk ratio (RRR, eb) says how many times higher the odds of one category rather than the baseline are per unit higher on a predictor, holding the others constant. Predicted probabilities turn the whole model back into chances. And a pseudo-R² summarizes how much the model improves on knowing nothing.

\text{RRR} = e^{b} = e^{0.0612} = 1.063, \qquad e^{10 \times 0.0612} = 1.84

Among students of the same gender, those one point higher in trait anxiety have 1.063 times the odds of avoiding rather than problem-solving, and those 10 points higher (about one SD) 1.84 times. At the same anxiety level, women's odds of seeking support rather than problem-solving are 2.19 times men's; unadjusted, the gender-by-strategy table gives 2.34, because women here score a little higher on anxiety. Both fall in the benchmarks' small range, but those cutoffs assume an outcome in about 1% of the comparison group, and every strategy here is common, so the labels understate them; the predicted probabilities are the better guide.

An RRR depends on the predictor's units, so a per-point value near 1 can still add up to a large change over the predictor's range. Judge a continuous predictor over a meaningful change, and always pair RRRs with predicted probabilities: 1.84 per 10 points compounds to about 3.4 over 20 points (1.84 × 1.84), and a man's predicted chance of avoiding is 16% at anxiety 35 against 35% at 55. The direction can even flip: anxiety's coefficient for seeking support is positive, yet women's predicted chance of seeking support falls from .44 at anxiety 45 to .38 at 74, the highest score in the sample, because avoiding gains faster still. McFadden's pseudo-R² of .039 says the model as a whole improves the log-likelihood only modestly; it has no agreed benchmarks and isn't a share of variance explained.

Benchmarks are conventions. Whether a difference of 16% against 35% in avoidant coping matters is a question about students and exams, not statistics.

Benchmarks Relative Risk Ratio (RRR) at a Glance
Below Small
0.60 < OR < 1.68
Small
1.68 ≤ OR < 3.47 (0.29 < OR ≤ 0.60)
Medium
3.47 ≤ OR < 6.71 (0.15 < OR ≤ 0.29)
Large
OR ≥ 6.71 (OR ≤ 0.15)

Full benchmarks and caveats ↑

Statistical Significance

A multinomial model answers three kinds of significance question. The model's likelihood-ratio test asks: if none of the predictors were related to the outcome, how often would adding them cut the deviance by this much or more? Each predictor's likelihood-ratio test asks the same about that predictor across all the comparisons at once. And each coefficient's z test asks whether one predictor is related to the odds of one category rather than the baseline. Below the conventional α = .05, a result is called statistically significant.

If anxiety and gender predicted nothing, a deviance drop of 33.41 or more would turn up less than once in 1,000 studies (p < .001). Anxiety across both comparisons: p < .001; gender: p = .005. Within the comparisons, anxiety's coefficient for avoiding has p < .001, but women's for avoiding has p = .306.

Test a predictor with its likelihood-ratio test first; its single z tests depend on which category is the baseline, and switching the baseline changes them. p is not the probability that a predictor is useless, and it says nothing about how much it matters. That's what the RRRs, their intervals and the predicted probabilities are for.

Benchmarks p at a Glance
Very Strong
p < .001
Strong
.001 ≤ p < .01
Significant
.01 ≤ p < .05
Weak
.05 ≤ p < .10
Little
p ≥ .10

Strength of evidence against the null hypothesis.

Full benchmarks and caveats ↑

Power and Sample Size

Power is the chance that a study detects an effect that's really there. A real effect doesn't guarantee a significant result: with too few people, it can easily come out non-significant. For a multinomial model, power depends on how strongly each predictor is related to the odds, how many people fall in each category (the smallest category matters most), how many predictors and categories there are, and α.

There's no simple formula or standard calculator for a multinomial model, so the usual approach is simulation: assume the effects you expect, generate many data sets of a given size from them, fit the model to each, and count how often the test you care about comes out significant.

Suppose the true effects were exactly the ones we found. Simulating 1,000 studies at each sample size, with students' anxiety and gender resampled from our data and their strategies drawn from the fitted model, the 2-df likelihood-ratio test for anxiety came out significant in about 70% of studies with 150 students, 84% with 200 and 99% with our 400; 80% power takes roughly 190. Gender's smaller effect needs more: about 75% with 300 students and 83% with 400, so roughly 360 for 80%. With 1,000 simulated studies, each percentage is within about 3 points of what endless simulations would give.

A single coefficient's 1-df test can have more or less power than its predictor's test across all comparisons: more when the effect sits mostly in one comparison, less when it's spread out. Here anxiety's avoid coefficient (p = .00001) and gender's support coefficient (p = .001) both beat their predictors' joint tests (p = .00003 and .005), so simulate the test you'll actually interpret. A rare category makes every comparison with it imprecise. Sample size also matters beyond significance: de Jong and colleagues (2019) found that multinomial models fitted by maximum likelihood are overfitted in small to medium samples, so their predictions look better in the sample than they will in new data; they stress having enough events (people in the smaller categories) per predictor, and recommend penalized estimation when building prediction models from small data sets.

Researchers usually aim for at least 80% power. Plan it before collecting data, by simulation, for the smallest effect that would matter and for the specific test you'll interpret, and make sure the rarest category will still have plenty of people.

Other Info

Follow-Up and Post Hoc Tests

Post Hoc Tests

  • Refitting with Another Baseline: Use it when you also want the comparison between the two non-baseline strategies (seeking support against avoiding). Its estimate is just the difference of the printed coefficients (anxiety: 0.0612 − 0.0288 = 0.032, RRR 1.033). Refit with one of them as the baseline to get its standard error and test (here z = 2.36, p = .018); the model and its fit don't change, only which pairs are printed.

References

  • Theil, H. (1969). A multinomial extension of the linear logit model. International Economic Review, 10(3), 251–259. https://doi.org/10.2307/2525642
  • Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied logistic regression (3rd ed.). Wiley. https://doi.org/10.1002/9781118548387
  • StataCorp. (2025). mlogit — Multinomial (polytomous) logistic regression. In Stata base reference manual (Release 19). Stata Press. https://www.stata.com/manuals/rmlogit.pdf
    “Thus, the exponentiated value of a coefficient is the relative-risk ratio for a one-unit change in the corresponding variable” (Remarks and examples)
    “The default is to choose the most frequent outcome.” (Options)
  • Hausman, J., & McFadden, D. (1984). Specification tests for the multinomial logit model. Econometrica, 52(5), 1219–1240. https://doi.org/10.2307/1910997
  • Cheng, S., & Long, J. S. (2007). Testing for IIA in the multinomial logit model. Sociological Methods & Research, 35(4), 583–600. https://doi.org/10.1177/0049124106292361
    “The authors conclude that tests of the IIA assumption that are based on the estimation of a restricted choice set are unsatisfactory for applied work.” (Abstract)
    “even in well-specified models, IIA tests often reject the assumption when the alternatives seem distinct and often fail to reject IIA when the alternatives can reasonably be viewed as close substitutes” (Abstract)
  • de Jong, V. M. T., Eijkemans, M. J. C., van Calster, B., Timmerman, D., Moons, K. G. M., Steyerberg, E. W., & van Smeden, M. (2019). Sample size considerations and predictive performance of multinomial logistic prediction models. Statistics in Medicine, 38(9), 1601–1619. https://doi.org/10.1002/sim.8063
    “It is shown that MLR estimated by Maximum Likelihood yields overfitted prediction models in small to medium sized data.” (Abstract)
  • Chen, H., Cohen, P., & Chen, S. (2010). How big is a big odds ratio? Interpreting the magnitudes of odds ratios in epidemiological studies. Communications in Statistics—Simulation and Computation, 39(4), 860–864. https://doi.org/10.1080/03610911003650383
    “OR = 1.68, 3.47, and 6.71 are equivalent to Cohen's d = 0.2 (small), 0.5 (medium), and 0.8 (large), respectively, when disease rate is 1% in the nonexposed group” (Abstract)
  • Wilks, S. S. (1938). The large-sample distribution of the likelihood ratio for testing composite hypotheses. The Annals of Mathematical Statistics, 9(1), 60–62. https://doi.org/10.1214/aoms/1177732360
  • McFadden, D. (1974). Conditional logit analysis of qualitative choice behavior. In P. Zarembka (Ed.), Frontiers in econometrics (pp. 105–142). Academic Press. https://eml.berkeley.edu/reprints/mcfadden/zarembka.pdf
    “we can define a coefficient of determination that is analogous to the multiple-correlation coefficient in the linear statistical model” (p. 121)
  • Nagelkerke, N. J. D. (1991). A note on a general definition of the coefficient of determination. Biometrika, 78(3), 691–692. https://doi.org/10.1093/biomet/78.3.691
    “A generalization of the coefficient of determination R2 to general regression models is discussed.” (Abstract)
  • Tjur, T. (2009). Coefficients of determination in logistic regression models—A new proposal: The coefficient of discrimination. The American Statistician, 63(4), 366–372. https://doi.org/10.1198/tast.2009.08210
    “namely the difference between the averages of fitted values for successes and failures, respectively” (Abstract)
  • Muff, S., Nilsen, E. B., O'Hara, R. B., & Nater, C. R. (2022). Rewriting results sections in the language of evidence. Trends in Ecology & Evolution, 37(3), 203–210. https://doi.org/10.1016/j.tree.2021.10.009
    “we suggest language of evidence that allows for a more nuanced approach to communicate scientific findings as a simple and intuitive alternative to statistical significance testing” (Abstract)
  • Box, G. E. P., & Tidwell, P. W. (1962). Transformation of the independent variables. Technometrics, 4(4), 531–550. https://doi.org/10.1080/00401706.1962.10490038
  • O'Brien, R. M. (2007). A caution regarding rules of thumb for variance inflation factors. Quality & Quantity, 41(5), 673–690. https://doi.org/10.1007/s11135-006-9018-6

Spot an error? Something still confusing? Tell us at support@statsforpsychology.com.