Markdown Editor

Logistic Regression Method Card Template (Statistics, R and Python)

Binary logistic regression: logit link, odds ratios, Wald vs likelihood-ratio tests, profile vs Wald CIs, predicted probabilities and marginal effects, classification metrics, calibration, separation and EPV, with R and Python

АкадемічніПовнийCox (1958); Nelder and Wedderburn (1972); Albert and Anderson (1984); Firth (1993); Hanley and McNeil (1982); Peduzzi et al. (1996); APA JARS-Quant (name)statslogisticregressionmath

Як використовувати: Use when the outcome is binary and you want adjusted odds ratios or predicted probabilities. Count events per variable first, check for separation, report odds ratios with profile CIs and likelihood-ratio tests, and translate at least one result to the probability scale. The main pitfall is reading an odds ratio as a risk ratio or a coefficient as a change in probability.

Перегляд

Logistic Regression: [Your study question or course]

A method card for binary logistic regression: when to use it, the logit link, odds and odds ratios, how to read each coefficient, Wald and likelihood-ratio tests, profile and Wald confidence intervals, predicted probabilities and average marginal effects, classification metrics with their threshold caveats, calibration, separation and events per variable. It is for students preparing for a course or exam, thesis writers with a yes/no outcome and analysts who build or check simple risk models. It speeds up learning by fixing the order of decisions, giving a fully worked example in R and Python, then fading the support in a half-solved example and an exercise, and testing you before you see the answers. The worked examples use fictional data; replace them with your own. Based on generalised linear model theory (Nelder and Wedderburn 1972; Cox 1958), the NIST/SEMATECH e-Handbook of Statistical Methods (public domain) and the method papers listed at the end.

Document control

FieldValue
Course or project[STAT 310 Generalised linear models / thesis chapter 5]
Learner or analyst[Your name]
Data set[Example: Lakemoor University pass study, n = 30, fictional]
Software and versions[Example: R 4.5.1 (stats); Python 3.14 with numpy 2.5.4, pandas 3.1.0, scipy 1.18.1, statsmodels 0.15.0]
DateYYYY-MM-DD
Version[0.1]

When to use it

One-sentence rule: use logistic regression when the outcome has exactly two categories (pass or fail, defect or no defect) and you want to model how the probability of the event changes with one or more predictors.

Use whenDo not use when
The outcome is binary and coded 0/1 (1 = the event)The outcome is continuous: use linear regression
You want odds ratios adjusted for other predictorsThe outcome is a count of events or a rate: use Poisson or negative binomial regression
You want predicted probabilities for new casesThe outcome has three or more categories: use multinomial or ordinal logistic regression
Observations are independentRows are clustered (students within classes, repeated tests per person): use a mixed or GEE logistic model
You have enough events for the number of predictorsOne predictor perfectly separates the two outcomes: use Firth's penalised logistic regression or exact methods

Alternatives: a chi-square or Fisher test for one categorical predictor without adjustment; the log-binomial or modified Poisson model when you want risk ratios rather than odds ratios (common outcomes, where odds ratios exaggerate); probit regression gives nearly the same fitted probabilities with a different link; classification trees or penalised models for pure prediction with many predictors.

The question and the data

Example question: at Lakemoor University (fictional), how do weekly study hours and attending the optional pre-exam review session relate to the probability of passing the introductory statistics exam?

VariableRoleMeasurement levelCoding
passedOutcomeBinary1 = passed, 0 = failed
hoursPredictorContinuous (ratio)Mean weekly study hours
reviewPredictorBinary1 = attended the review session, 0 = did not (reference)

Design: cross-sectional observational study, sample size n = 30 students, 19 passed and 11 failed. The smaller of the two counts (11) is the number that limits model size: with 2 predictors that is 11 / 2 = 5.5 events per variable (EPV), below the traditional rule of thumb of 10 (Peduzzi et al. 1996). The example is fine for learning, but in a real study it would be a reason for caution about the precision of the estimates.

StudenthoursreviewpassedStudenthoursreviewpassed
110016710
220017801
321018811
430019811
531120900
640021911
7400221001
8510231011
9501241101
10511251111
11600261211
12611271201
13600281311
14711291401
15701301511

The chart shows the outcome (0 or 1, slightly jittered by colour) against study hours, with the fitted probability of passing from the one-predictor model passed ~ hours (Example B). The S-shaped curve is the point of the logit link: probabilities stay between 0 and 1, and the effect of one more hour on the probability is largest near p = 0.5.

Assumptions and how to check each

AssumptionHow to checkRPythonIf violated
Binary outcome, one row per independent unitStudy design; table(y)table(d$passed)d["passed"].value_counts()Clustered data: mixed logistic (lme4::glmer) or GEE
Linearity of the log-odds in each continuous predictorAdd a squared term or spline and compare by likelihood ratio; binned residual plotanova(fit, update(fit, . ~ . + I(hours^2)), test = "LRT")Compare llf of the two fits; 2 * (llf1 - llf0)Spline or polynomial term, or a transformation (log of a skewed predictor)
Enough events per variableCount the smaller outcome group and divide by the number of coefficients (not counting the intercept)min(table(d$passed)) / 2d["passed"].value_counts().min() / 2Fewer predictors, penalised (Firth) regression, collect more data; for prediction models use formal sample size methods
No (quasi-)complete separationHuge coefficients and SEs, warning "fitted probabilities numerically 0 or 1"Read warnings; summary(fit)PerfectSeparationWarning or huge bseFirth regression (logistf package), exact logistic, or merge sparse categories
No severe multicollinearityVIF on the same predictorscar::vif(fit)variance_inflation_factor on the design matrixAs for linear regression
No overly influential pointsCook's distance and dfbeta for GLMsplot(fit, which = 5); influence.measures(fit)fit_glm.get_influence().summary_frame()Check for errors, report with and without

Notice what is not in the table: residuals are not expected to be normal, and their variance is not constant (it is p(1−p)p(1-p) by design). Normality tests on logistic residuals are a common wasted step.

How it works

In plain language. A straight line for a probability would give values below 0 and above 1. Logistic regression instead models the log-odds of the event as a straight line, and then turns the log-odds back into a probability with an S-shaped curve. Each coefficient is the change in log-odds per unit of the predictor, holding the others constant. Exponentiating it gives an odds ratio: the factor by which the odds of the event are multiplied for a one-unit increase. The coefficients are found by maximum likelihood: the values that make the observed passes and fails most probable.

Three scales. For probability pp:

odds=p1−p,logit⁡(p)=ln⁡ ⁣(p1−p),p=11+e−η\text{odds} = \frac{p}{1 - p}, \qquad \operatorname{logit}(p) = \ln\!\left(\frac{p}{1 - p}\right), \qquad p = \frac{1}{1 + e^{-\eta}}

Model for observation ii, with linear predictor ηi\eta_i:

ln⁡ ⁣(pi1−pi)=ηi=β0+β1xi1+⋯+βkxik,Yi∼Bernoulli(pi)\ln\!\left(\frac{p_i}{1 - p_i}\right) = \eta_i = \beta_0 + \beta_1 x_{i1} + \cdots + \beta_k x_{ik}, \qquad Y_i \sim \text{Bernoulli}(p_i)

Odds ratio for a one-unit increase in xjx_j, and for a cc-unit increase:

ORj=eβj,ORj(c)=ec βj\text{OR}_j = e^{\beta_j}, \qquad \text{OR}_j(c) = e^{c\,\beta_j}

Log-likelihood maximised by the software (iteratively reweighted least squares):

ℓ(β)=∑i=1n[yiln⁡pi+(1−yi)ln⁡(1−pi)]\ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \left[ y_i \ln p_i + (1 - y_i)\ln(1 - p_i) \right]

Tests. Wald test of one coefficient, and likelihood-ratio (LR) test of a full model against a reduced one with qq fewer coefficients, where DD is the deviance −2ℓ-2\ell:

z=β^jSE(β^j),G2=Dreduced−Dfull=2(ℓfull−ℓreduced)∼χq2z = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)}, \qquad G^2 = D_{\text{reduced}} - D_{\text{full}} = 2(\ell_{\text{full}} - \ell_{\text{reduced}}) \sim \chi^2_q

Average marginal effect (AME) of a continuous predictor on the probability scale, and of a binary predictor (discrete change):

AMEj=1n∑i=1nβ^j p^i(1−p^i),AMEbinary=1n∑i=1n[p^i(x=1)−p^i(x=0)]\text{AME}_j = \frac{1}{n}\sum_{i=1}^{n} \hat{\beta}_j\, \hat{p}_i (1 - \hat{p}_i), \qquad \text{AME}_{\text{binary}} = \frac{1}{n}\sum_{i=1}^{n} \left[\hat{p}_i(x = 1) - \hat{p}_i(x = 0)\right]
SymbolMeaning
pip_iProbability of the event (Y = 1) for observation ii
ηi\eta_iLinear predictor, the log-odds of the event
βj\beta_jChange in log-odds per unit of xjx_j, others held constant
eβje^{\beta_j}Odds ratio per unit of xjx_j
ℓ\ell, DDLog-likelihood and deviance (D=−2ℓD = -2\ell for binary data)
G2G^2Likelihood-ratio statistic, compared with χ2\chi^2 on qq df
AMEAverage change in predicted probability, in probability units
EPVEvents per variable: smaller outcome count divided by the number of coefficients

Profile vs Wald intervals. The Wald interval is β^±1.96 SE\hat{\beta} \pm 1.96\,SE, symmetric on the log-odds scale. The profile likelihood interval follows the actual shape of the likelihood and is usually more accurate in small samples or near separation. In R, confint() on a glm gives profile intervals (since R 4.4.0 the method lives in stats; earlier versions borrowed it from MASS), while confint.default() gives Wald intervals. In statsmodels, conf_int() gives Wald intervals. So R and Python will print different CIs for the same model, and both are "correct" for their method.

Classification. Turning probabilities into yes/no predictions needs a threshold tt: predict 1 when p^≥t\hat{p} \ge t. Sensitivity (true positive rate) and specificity (true negative rate) depend on tt; the ROC curve plots them for every threshold, and the area under it (AUC, or the c-statistic) is the probability that a randomly chosen event has a higher predicted probability than a randomly chosen non-event. Calibration asks a different question: among cases given about 70%, do about 70% have the event?

Worked example A (full): passing the exam

Hypotheses. For each predictor, H0:βj=0H_0: \beta_j = 0 (odds ratio 1, no association with passing once the other predictor is held constant) against H1:βj≠0H_1: \beta_j \neq 0, two-sided, α = 0.05. Overall: H0H_0: neither predictor is related to passing.

Checks before reading the results. Events: 19 passes and 11 fails, so EPV = 5.5 (low; noted as a limitation). No separation warning appears, and the standard errors are of normal size. Adding I(hours^2) lowers the deviance by only 0.06 (LR test p = .803), so log-odds linearity in hours is acceptable.

R code.

d <- data.frame(
  hours  = c(1, 2, 2, 3, 3, 4, 4, 5, 5, 5, 6, 6, 6, 7, 7,
             7, 8, 8, 8, 9, 9, 10, 10, 11, 11, 12, 12, 13, 14, 15),
  review = c(0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0,
             1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1),
  passed = c(0, 0, 0, 0, 1, 0, 0, 0, 1, 1, 0, 1, 0, 1, 1,
             0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1))

fit <- glm(passed ~ hours + review, family = binomial, data = d)
summary(fit)
confint(fit)                       # profile likelihood CIs (stats::confint, R >= 4.4)
confint.default(fit)               # Wald CIs: estimate +/- 1.96 SE
round(exp(cbind(OR = coef(fit), confint(fit))), 3)
drop1(fit, test = "LRT")           # likelihood-ratio test for each term

# Predicted probabilities and average marginal effects (AME)
p <- fitted(fit)
predict(fit, data.frame(hours = c(4, 8), review = 1), type = "response")
mean(p * (1 - p)) * coef(fit)["hours"]                    # AME of one more hour
d1 <- transform(d, review = 1); d0 <- transform(d, review = 0)
mean(predict(fit, d1, type = "response") - predict(fit, d0, type = "response"))

# Classification at threshold 0.5, and AUC via the rank (Mann-Whitney) formula
table(predicted = as.integer(p >= 0.5), actual = d$passed)
r <- rank(p); n1 <- sum(d$passed); n0 <- sum(1 - d$passed)
(sum(r[d$passed == 1]) - n1 * (n1 + 1) / 2) / (n1 * n0)

# Calibration: observed vs mean predicted in three risk groups
grp <- cut(p, c(0, 1/3, 2/3, 1))
aggregate(cbind(observed = d$passed, predicted = p) ~ grp, FUN = mean)

R output (R 4.5.1; output abridged: some calls are shown with a > prompt and a short comment so you can see which line printed what, and R's "Waiting for profiling to be done..." messages are omitted).

Coefficients:
            Estimate Std. Error z value Pr(>|z|)   
(Intercept)  -4.1784     1.6771  -2.491  0.01272 * 
hours         0.5961     0.2273   2.623  0.00871 **
review        1.7397     1.1060   1.573  0.11572   

    Null deviance: 39.429  on 29  degrees of freedom
Residual deviance: 22.248  on 27  degrees of freedom
AIC: 28.248

> confint(fit)            # profile
                 2.5 %    97.5 %
(Intercept) -8.2697881 -1.441756
hours        0.2326828  1.159566
review      -0.2981269  4.211865
> confint.default(fit)    # Wald
                 2.5 %     97.5 %
(Intercept) -7.4654175 -0.8913465
hours        0.1507225  1.0415589
review      -0.4280001  3.9073975

               OR 2.5 % 97.5 %
(Intercept) 0.015 0.000  0.237
hours       1.815 1.262  3.189
review      5.696 0.742 67.482

       Df Deviance    AIC     LRT  Pr(>Chi)    
<none>      22.248 28.248                      
hours   1   35.740 39.740 13.4918 0.0002396 ***
review  1   25.021 29.021  2.7735 0.0958393 .  

> predict(fit, data.frame(hours = c(4, 8), review = 1), type = "response")
        1         2 
0.4864732 0.9113671 
> mean(p * (1 - p)) * coef(fit)["hours"]
     hours 
0.06996238 
> mean(predict(fit, d1, ...) - predict(fit, d0, ...))
[1] 0.2263289

         actual
predicted  0  1
        0  8  3
        1  3 16
[1] 0.9043062

            grp  observed predicted
1     (0,0.333] 0.1428571 0.1284975
2 (0.333,0.667] 0.5714286 0.4939522
3     (0.667,1] 0.8750000 0.9151782

Python code.

import numpy as np
import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf
from scipy.stats import mannwhitneyu

d = pd.DataFrame({
    "hours":  [1, 2, 2, 3, 3, 4, 4, 5, 5, 5, 6, 6, 6, 7, 7,
               7, 8, 8, 8, 9, 9, 10, 10, 11, 11, 12, 12, 13, 14, 15],
    "review": [0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0,
               1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1],
    "passed": [0, 0, 0, 0, 1, 0, 0, 0, 1, 1, 0, 1, 0, 1, 1,
               0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]})

fit = smf.logit("passed ~ hours + review", data=d).fit(disp=0)
print(fit.summary().tables[1])
print("LR chi2 =", round(fit.llr, 3), " p =", round(fit.llr_pvalue, 7),
      " McFadden R2 =", round(fit.prsquared, 4))
print("OR:", np.exp(fit.params).round(3).to_dict())
print("Wald 95% CI for OR:\n", np.exp(fit.conf_int()).round(3))

# Average marginal effects; dummy=True gives the discrete 0 -> 1 change for review
print(fit.get_margeff(at="overall", dummy=True).summary().tables[1])

# Same model as a GLM (identical estimates; deviance and AIC as in R)
glm = smf.glm("passed ~ hours + review", data=d, family=sm.families.Binomial()).fit()
print("GLM deviance =", round(glm.deviance, 3), " AIC =", round(glm.aic, 3))

p = fit.predict(d)
print(pd.crosstab((p >= 0.5).astype(int), d["passed"], rownames=["predicted"]))
u = mannwhitneyu(p[d.passed == 1], p[d.passed == 0]).statistic
print("AUC =", round(u / (19 * 11), 4))

Python output (Python 3.14, statsmodels 0.15.0, scipy 1.18.1; output abridged: the ==== table borders are omitted).

                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept     -4.1784      1.677     -2.491      0.013      -7.465      -0.891
hours          0.5961      0.227      2.623      0.009       0.151       1.042
review         1.7397      1.106      1.573      0.116      -0.428       3.907
LR chi2 = 17.181  p = 0.0001858  McFadden R2 = 0.4358
OR: {'Intercept': 0.015, 'hours': 1.815, 'review': 5.696}
Wald 95% CI for OR:
                0       1
Intercept  0.001   0.410
hours      1.163   2.834
review     0.652  49.769
                dy/dx    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
hours          0.0700      0.011      6.176      0.000       0.048       0.092
review         0.2263      0.136      1.662      0.096      -0.041       0.493
GLM deviance = 22.248  AIC = 28.248
passed     0   1
predicted
0          8   3
1          3  16
AUC = 0.9043

Where R and Python differ, and why. The estimates, SEs, Wald z tests, deviance and AIC agree exactly. The confidence intervals do not: R's confint() gives profile intervals (OR for review [0.742, 67.482]) while statsmodels gives Wald intervals (OR [0.652, 49.769]); the Wald intervals in R (confint.default) match Python. The marginal effect of review also depends on a setting: with dummy=True statsmodels computes the discrete change from 0 to 1 (0.2263, the same as the R calculation); without it, statsmodels treats review as continuous and reports a derivative of 0.2042. For a 0/1 predictor the discrete change is the one to report.

Calculation by hand for the key statistics.

  • Odds ratio for hours: e0.5961=1.815e^{0.5961} = 1.815. Each extra weekly hour multiplies the odds of passing by 1.82, an 82% increase in the odds.
  • For 3 extra hours: e3×0.5961=e1.788=5.98e^{3 \times 0.5961} = e^{1.788} = 5.98.
  • Predicted probability for a review attender who studies 4 hours: η=−4.1784+0.5961×4+1.7397=−0.0543\eta = -4.1784 + 0.5961 \times 4 + 1.7397 = -0.0543, so p=1/(1+e0.0543)=0.486p = 1/(1 + e^{0.0543}) = 0.486.
  • LR test for the whole model: G2=39.429−22.248=17.181G^2 = 39.429 - 22.248 = 17.181 on 2 df, p<.001p < .001.
  • From the confusion matrix (threshold 0.5): sensitivity =16/19=0.842= 16/19 = 0.842, specificity =8/11=0.727= 8/11 = 0.727, accuracy =24/30=0.80= 24/30 = 0.80.

Wald vs likelihood-ratio for review. The Wald test gives p = .116, the LR test p = .096. They disagree because the Wald test relies on a symmetric normal approximation that is poor when a coefficient is large relative to the information in the data (here only 30 students and a wide, skewed likelihood). In small samples prefer the LR test and the profile CI; both say the review effect is imprecisely estimated (OR 5.70, profile 95% CI [0.74, 67.48]).

Effect size with CI. Odds ratios with profile CIs: hours OR = 1.82, 95% CI [1.26, 3.19]; review OR = 5.70, 95% CI [0.74, 67.48]. On the probability scale (AMEs), one more weekly hour raises the probability of passing by about 7.0 percentage points on average (95% CI [4.8, 9.2] from the delta method in statsmodels), and attending the review session by about 22.6 points (95% CI [-4.1, 49.3]). Discrimination: AUC = 0.90.

Interpretation. Holding review attendance constant, students who study more are much more likely to pass: the odds of passing rise by a factor of 1.82 for each extra weekly hour. Review attendance is associated with higher odds of passing, but the interval is very wide and includes 1, so these data cannot tell a large benefit from none. The model separates passers from failers well (AUC 0.90) and is roughly calibrated in three risk groups (observed 14%, 57% and 88% against predicted 13%, 49% and 92%), but with EPV = 5.5 and no external validation, it should not be used to predict outcomes for other cohorts.

APA-style report. A logistic regression predicted passing the exam from weekly study hours and review session attendance, χ²(2) = 17.18, p < .001, McFadden R² = .44. Each additional study hour was associated with higher odds of passing, OR = 1.82, 95% CI [1.26, 3.19] (profile), LR χ²(1) = 13.49, p < .001; the average marginal effect was 7.0 percentage points per hour. Review attendance was not clearly associated with passing, OR = 5.70, 95% CI [0.74, 67.48], LR χ²(1) = 2.77, p = .096. The model discriminated well (AUC = .90); with 11 failures and 2 predictors, estimates are imprecise.

Worked example B (half solved): one predictor

Model passed ~ hours on the same 30 students. R output (abridged: coef(summary(fit1)), and the last line summarises exp(confint(fit1)) in words):

              Estimate Std. Error   z value    Pr(>|z|)
(Intercept) -3.2632042  1.3955360 -2.338316 0.019370860
hours        0.5798946  0.2155692  2.690062 0.007143865
profile 95% CI for the odds ratio of hours: [1.264, 3.039]
  1. Odds ratio for hours: e0.5799=e^{0.5799} = [ ].
  2. In words: each extra weekly hour multiplies the odds of passing by [ ], that is, the odds increase by [ ]%.
  3. Log-odds for a student with 5 hours: η=−3.2632+0.5799×5=\eta = -3.2632 + 0.5799 \times 5 = [ ].
  4. Probability: p=1/(1+e−η)=p = 1/(1 + e^{-\eta}) = [ ].
  5. The hours at which the predicted probability is exactly 0.5 (where η=0\eta = 0): 3.2632/0.5799=3.2632 / 0.5799 = [ ] hours.
  6. Report sentence: OR = [ ], 95% CI [ ], z = [ ], p = [ ].

Exercise C (unsolved)

Brightwell Bakery (fictional) records, for 20 batches of rye loaves, how far the oven temperature was from the set point (°C, negative = too cold) and whether the batch failed the quality check (defect = 1). The baker suspects that deviations in either direction cause defects.

Batchtemp_devdefectBatchtemp_devdefect
1-1211100
2-1011210
3-901321
4-711430
5-601540
6-501650
7-411761
8-301881
9-2019101
10-1020121

Tasks: (a) explain why defect ~ temp_dev would be the wrong model for the baker's idea, and create abs_dev = abs(temp_dev); (b) fit defect ~ abs_dev; report the odds ratio per °C with its profile 95% CI and the LR test; (c) give the predicted probability of a defect at 0 °C and at 10 °C of deviation; (d) compute EPV and comment.

Reading software output

Output line (R summary(glm))statsmodels equivalentMeaningReport it?
EstimatecoefChange in log-odds per unit, others held constantUsually as OR = exp(estimate)
Std. Error, z value, Pr(>|z|)std err, z, P>|z|Wald test of βj=0\beta_j = 0Prefer the LR test in small samples
confint(fit)not built in (Wald only)Profile likelihood CI for βj\beta_j; exponentiate for ORYes, say "profile"
confint.default(fit)conf_int()Wald CISay "Wald" if you use it
Null deviance, Residual deviance-2 * llnull, -2 * llf (or GLM deviance)Model fit; their difference is the LR statisticYes, as χ²(df)
AICaicComparison of non-nested models on the same dataWhen comparing models
drop1(fit, test = "LRT")compare llf of nested fitsLR test for each termYes
no direct lineprsquaredMcFadden pseudo-R², not comparable with OLS R²Optional, name the type
Number of Fisher Scoring iterationsconvergedConvergence; very many iterations can signal separationCheck, do not report

Reporting template

Fill-in sentence. A logistic regression modelled the probability of [event] from [predictors] (n = [n], [events] events), χ²([df]) = [G²], p [=/<] [p]. Holding [other predictors] constant, each [unit] increase in [predictor] was associated with [higher/lower] odds of [event], OR = [OR], 95% CI [lower, upper] ([profile/Wald]), LR χ²(1) = [value], p [=/<] [p]; the average marginal effect was [AME] percentage points. [For prediction: AUC = [value]; calibration [result]; threshold [t] gave sensitivity [ ] and specificity [ ].]

PredictorbSEOR95% CI for ORLR χ²(1)p
Example: Study hours (per hour)0.5960.2271.82[1.26, 3.19]13.49< .001
Example: Review attended (vs not)1.7401.1065.70[0.74, 67.48]2.77.096
[Predictor][ ][ ][ ][ ][ ][ ]

Results checklist.

  • Event coded 1 and reference levels stated
  • Number of events and non-events, and EPV, reported
  • Separation checked (no warnings, no huge SEs)
  • Odds ratios with CIs, CI method named (profile or Wald)
  • LR tests used, or Wald tests justified
  • At least one result translated to probabilities (predicted probabilities or AMEs)
  • For prediction aims: AUC, calibration, threshold and validation reported
  • Odds ratios not described as risk ratios
  • Software and versions named

Common mistakes and misconceptions

  • "OR = 1.82 means passing is 1.82 times as likely." That is a risk ratio. An odds ratio multiplies the odds. When the event is common (here 63% passed), the odds ratio is further from 1 than the risk ratio.
  • Reading a coefficient as a change in probability. 0.596 is a change in log-odds. The change in probability depends on where you start: it is largest near p = 0.5 and small near 0 or 1. Report AMEs or predicted probabilities for that.
  • Comparing odds ratios or coefficients across models with different predictors as if only the predictor of interest changed. In logistic regression adding a predictor that is unrelated to the others (but related to the outcome) still changes the other coefficients (non-collapsibility).
  • Trusting Wald p-values when SEs are large. Use LR tests and profile intervals, especially in small samples.
  • Ignoring separation warnings. In a small perfectly separated example (8 students, all with 5 or more hours passed), R reports an estimate of 45.8 with SE 80,644 and p = 1.00. That is not "no effect": the maximum likelihood estimate does not exist. Use Firth regression.
  • Reporting accuracy alone. With 63% passers, predicting "pass" for everyone already gives 63% accuracy. Report sensitivity and specificity, and say which threshold you used and why.
  • Treating the 0.5 threshold as part of the model. The threshold is a decision about costs, not a statistical result. A tutor who wants to catch every student at risk will choose a different threshold from one who has few tutoring places.
  • Using a high AUC as evidence of good calibration. AUC measures ranking only. A model can rank perfectly and still predict 90% for a group where 60% pass.
  • Too many predictors for the events. With 11 failures, a model with 8 predictors is badly overfitted, whatever its apparent AUC.

Self-test questions

  1. Convert a probability of 0.80 to odds and to log-odds.
  2. The coefficient of hours is 0.5961. State its meaning on the log-odds scale and on the odds scale in one sentence each.
  3. Why do R's confint(fit) and statsmodels' conf_int() give different intervals for the same logistic model, and which would you report for n = 30?
  4. With a threshold of 0.7 instead of 0.5, the confusion matrix becomes: predicted 0 / actual 0 = 9, predicted 0 / actual 1 = 5, predicted 1 / actual 0 = 2, predicted 1 / actual 1 = 14. Compute sensitivity and specificity, and say what changed.
  5. You fit a model and see the warning "fitted probabilities numerically 0 or 1 occurred" and a coefficient of 38 with SE 9,000. What is happening and what do you do?
  6. A study has 400 participants and 24 events. How many predictors can you reasonably include by the EPV = 10 guide, and what else should you consider?

My notes and next steps

ItemEntry
Confidence with this method (1-5)[3]
What I can now do[Fit a logistic model and explain an odds ratio]
Still unclear[When to use risk ratios instead]
Mistakes to log in academic/stats-error-log[Called an OR a risk ratio]
Next review date (spacing)YYYY-MM-DD
Next template[academic/stats-count-regression]

Sources and further reading

Structure aligned with the method card used across this statistics library; formulas are standard generalised linear model results in our own notation.

  • Cox, D. R. (1958). The regression analysis of binary sequences. Journal of the Royal Statistical Society, Series B, 20(2), 215-232.
  • Nelder, J. A., and Wedderburn, R. W. M. (1972). Generalized linear models. Journal of the Royal Statistical Society, Series A, 135(3), 370-384.
  • Albert, A., and Anderson, J. A. (1984). On the existence of maximum likelihood estimates in logistic regression models. Biometrika, 71(1), 1-10.
  • Firth, D. (1993). Bias reduction of maximum likelihood estimates. Biometrika, 80(1), 27-38.
  • Hanley, J. A., and McNeil, B. J. (1982). The meaning and use of the area under a receiver operating characteristic (ROC) curve. Radiology, 143(1), 29-36.
  • Peduzzi, P., Concato, J., Kemper, E., Holford, T. R., and Feinstein, A. R. (1996). A simulation study of the number of events per variable in logistic regression analysis. Journal of Clinical Epidemiology, 49(12), 1373-1379.
  • Van Calster, B., McLernon, D. J., van Smeden, M., Wynants, L., and Steyerberg, E. W. (2019). Calibration: the Achilles heel of predictive analytics. BMC Medicine, 17, 230.
  • Riley, R. D., Ensor, J., Snell, K. I. E., and colleagues (2020). Calculating the sample size required for developing a clinical prediction model. BMJ, 368, m441.
  • Hosmer, D. W., Lemeshow, S., and Sturdivant, R. X. (2013). Applied Logistic Regression (3rd ed.). Wiley.
  • NIST/SEMATECH e-Handbook of Statistical Methods, https://www.itl.nist.gov/div898/handbook/ (accessed 2026-10-25). Public domain.
  • Further reading by name: de Smith, M. J., Statistical Analysis Handbook (StatsRef, 2024 edition); OpenIntro Statistics; APA Publication Manual (7th ed.) and APA JARS-Quant; Wasserstein and Lazar (2016), The ASA statement on p-values.

Answer key

Example B.

  1. e0.5799=1.786e^{0.5799} = 1.786.
  2. By 1.79; the odds increase by about 79%.
  3. −3.2632+2.8995=−0.3637-3.2632 + 2.8995 = -0.3637.
  4. 1/(1+e0.3637)=0.4101/(1 + e^{0.3637}) = 0.410 (R predict(..., type = "response"): 0.4100566).
  5. 3.2632/0.5799=5.633.2632 / 0.5799 = 5.63 hours.
  6. OR = 1.79, 95% CI [1.26, 3.04] (profile), z = 2.69, p = .007.

Exercise C (R 4.5.1; output abridged and relabelled: the values come from coef(summary()), exp(cbind(coef(), confint())), anova(fit, test = "LRT") and predict(type = "response"), condensed to one line each):

              Estimate Std. Error   z value   Pr(>|z|)
(Intercept) -2.7269241  1.2355324 -2.207084 0.02730817
abs_dev      0.4542173  0.2005087  2.265325 0.02349277
odds ratio and profile 95% CI:  abs_dev 1.57494013  [1.134215917, 2.593523]
LR test: Deviance 8.2491 on 1 df, Resid. Dev 19.276, Pr(>Chi) 0.004077
predicted probability at abs_dev 0: 0.0614032   at abs_dev 10: 0.8599950

(a) temp_dev assumes the log-odds rise steadily from too cold to too hot, so cold and hot deviations would cancel out; the absolute deviation (or a squared term) matches a U-shaped risk. (b) OR = 1.57 per °C of deviation, 95% CI [1.13, 2.59], LR χ²(1) = 8.25, p = .004. (c) 6.1% at the set point, 86.0% at 10 °C off. (d) 9 defects and 11 good batches, so the smaller group is 9 and EPV = 9 / 1 = 9, close to the guide of 10; one predictor is the most these data support.

Self-test.

  1. Odds =0.80/0.20=4= 0.80 / 0.20 = 4; log-odds =ln⁡4=1.386= \ln 4 = 1.386.
  2. Log-odds: each extra weekly hour raises the log-odds of passing by 0.596, review attendance held constant. Odds: each extra hour multiplies the odds of passing by e0.596=1.82e^{0.596} = 1.82.
  3. R's confint on a glm gives profile likelihood intervals; statsmodels gives Wald intervals. For n = 30, report the profile interval (more accurate when the likelihood is not symmetric) and name the method.
  4. Sensitivity =14/19=0.737= 14/19 = 0.737, specificity =9/11=0.818= 9/11 = 0.818. Raising the threshold makes the model predict "pass" less often, so it misses more passers (sensitivity falls from 0.842) and wrongly predicts fewer failers as passers (specificity rises from 0.727).
  5. Complete or quasi-complete separation: some combination of predictors predicts the outcome perfectly, so the maximum likelihood estimate goes to infinity. Look for the predictor or category that separates, then use Firth penalised regression (logistf in R) or exact logistic regression, or merge sparse categories if that is substantively justified.
  6. 24 / 10 = about 2 predictors (counting dummies and interaction terms as separate coefficients). EPV is only a rough guide: for prediction models, use a formal sample size method (Riley and colleagues, 2020), and consider penalised estimation.

Revision history

VersionDateAuthorChange
0.1YYYY-MM-DD[Name]First draft from the template
[0.2]YYYY-MM-DD[Name][Switched to profile CIs and LR tests]

Схожі шаблони