ML Arc, Foundations — Generative vs Discriminative (R companion)¶
e1071 (Naive Bayes), MASS (linear discriminant analysis, LDA), and glm (logistic) on the credit data¶
R's classic classification packages line up with the two philosophies: e1071::naiveBayes and MASS::lda are generative (they model $p(x\mid y)p(y)$), glm(family=binomial) is discriminative (it models $p(y\mid x)$). This companion reproduces the Python notebook's conclusions: Naive Bayes ranks well (AUC) but is badly miscalibrated, LDA sits in between, and a learning curve shows the generative model converging faster at small $n$.
Its numbers differ from the Python side, and the reason is worth chasing rather than waving at, because the obvious suspect is wrong. It is not the implementations. Run e1071::naiveBayes on Python's exact train/test indices and it returns AUC 0.7190, accuracy 0.5191 and ECE 0.3785 against Python's 0.7190, 0.5194 and 0.3781 — agreement to three decimals. The difference is entirely the split, and the fact that a split can move these numbers so far is itself the notebook's point, as section 1 shows.
1. The three classifiers — ranking vs calibration¶
We fit Naive Bayes, LDA, and logistic regression and score each by AUC (ranking), accuracy, and a simple calibration error. As in Python, Naive Bayes' AUC is competitive but its accuracy is far lower and its calibration far worse — the independence assumption gives good discrimination and bad probabilities. LDA, with a shared covariance instead of full independence, is much better calibrated.
suppressMessages({library(e1071); library(MASS)})
d<-read.csv("credit_default.csv"); d$default<-factor(d$default)
set.seed(0); i<-sample(nrow(d), 0.7*nrow(d)); tr<-d[i,]; te<-d[-i,]
feat<-setdiff(names(d),"default")
auc<-function(p,y){ y<-as.integer(y)-1; r<-rank(p); (sum(r[y==1])-sum(y==1)*(sum(y==1)+1)/2)/(sum(y==1)*sum(y==0)) }
ece<-function(p,y){ y<-as.integer(y)-1; b<-cut(p,seq(0,1,0.1),include.lowest=TRUE); e<-0
for(l in levels(b)){ m<-b==l; if(sum(m)>0) e<-e+mean(m)*abs(mean(y[m])-mean(p[m])) }; e }
nb<-naiveBayes(default~., data=tr); pnb<-predict(nb, te, type="raw")[,2]
ld<-lda(default~., data=tr); pld<-predict(ld, te)$posterior[,2]
lg<-glm(default~., data=tr, family=binomial); plg<-predict(lg, te, type="response")
res<-rbind(
`Naive Bayes (generative)`=c(AUC=auc(pnb,te$default), accuracy=mean((pnb>=.5)==(te$default=="1")), ECE=ece(pnb,te$default)),
`LDA (generative)` =c(auc(pld,te$default), mean((pld>=.5)==(te$default=="1")), ece(pld,te$default)),
`Logistic (discriminative)`=c(auc(plg,te$default), mean((plg>=.5)==(te$default=="1")), ece(plg,te$default)))
print(round(res,3))
cat("\nNaive Bayes: competitive AUC but low accuracy and high ECE (miscalibrated); LDA better calibrated; logistic best accuracy.\n")
# How much does the SPLIT move these? Naive Bayes' accuracy turns out to be remarkably fragile.
set.seed(99); nr <- 20
sa <- sacc <- sece <- numeric(nr); la <- lacc <- numeric(nr)
for(s in 1:nr){
ii <- sample(nrow(d), 0.7*nrow(d)); t1 <- d[ii,]; t2 <- d[-ii,]
p2 <- predict(naiveBayes(default~., data=t1), t2, type="raw")[,2]
g2 <- suppressWarnings(glm(default~., data=t1, family=binomial))
q2 <- predict(g2, t2, type="response")
sa[s]<-auc(p2,t2$default); sacc[s]<-mean((p2>=.5)==(t2$default=="1")); sece[s]<-ece(p2,t2$default)
la[s]<-auc(q2,t2$default); lacc[s]<-mean((q2>=.5)==(t2$default=="1"))
}
cat(sprintf("\nOver %d random 70/30 splits:\n", nr))
cat(sprintf(" Naive Bayes AUC %.3f (sd %.3f) accuracy %.3f (sd %.3f) ECE %.3f (sd %.3f)\n",
mean(sa), sd(sa), mean(sacc), sd(sacc), mean(sece), sd(sece)))
cat(sprintf(" Logistic AUC %.3f (sd %.3f) accuracy %.3f (sd %.3f)\n", mean(la), sd(la), mean(lacc), sd(lacc)))
cat(sprintf("\nNaive Bayes' AUC barely moves (sd %.3f) -- its RANKING is stable. Its accuracy swings with sd %.3f,\n", sd(sa), sd(sacc)))
cat(sprintf("about %.0f times logistic's %.3f. That is not a second defect; it is the SAME defect seen from another angle.\n", sd(sacc)/sd(lacc), sd(lacc)))
cat("A model whose probabilities pile up at 0 and 1 has an accuracy that hinges on exactly how many observations\n")
cat("fall either side of the 0.5 cut, so a different split moves it a long way. Logistic spreads its probabilities\n")
cat("across the interval, so the same threshold lands in a sparse region and the accuracy is stable.\n")
cat("\nThis also explains why this notebook's table differs from the Python one: same data, same algorithm, different\n")
cat("split. Given Python's exact indices, e1071 reproduces its numbers to three decimals. The ranking metric agrees\n")
cat("across splits and the threshold metric does not, which is the whole lesson about miscalibration in one line.\n")
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
AUC accuracy ECE Naive Bayes (generative) 0.735 0.703 0.255 LDA (generative) 0.705 0.812 0.045 Logistic (discriminative) 0.712 0.810 0.049
Naive Bayes: competitive AUC but low accuracy and high ECE (miscalibrated); LDA better calibrated; logistic best accuracy.
Over 20 random 70/30 splits:
Naive Bayes AUC 0.736 (sd 0.007) accuracy 0.687 (sd 0.069) ECE 0.253 (sd 0.044)
Logistic AUC 0.722 (sd 0.004) accuracy 0.809 (sd 0.004)
Naive Bayes' AUC barely moves (sd 0.007) -- its RANKING is stable. Its accuracy swings with sd 0.069,
about 16 times logistic's 0.004. That is not a second defect; it is the SAME defect seen from another angle.
A model whose probabilities pile up at 0 and 1 has an accuracy that hinges on exactly how many observations
fall either side of the 0.5 cut, so a different split moves it a long way. Logistic spreads its probabilities
across the interval, so the same threshold lands in a sparse region and the accuracy is stable.
This also explains why this notebook's table differs from the Python one: same data, same algorithm, different
split. Given Python's exact indices, e1071 reproduces its numbers to three decimals. The ranking metric agrees
across splits and the threshold metric does not, which is the whole lesson about miscalibration in one line.
2. The Ng-Jordan learning curve — generative converges faster¶
On a controlled dataset where Naive Bayes' assumptions hold (independent Gaussian features), we plot test error against training size. Naive Bayes wins when data are scarce; logistic regression catches up as $n$ grows.
One caveat has to be attached to the small-$n$ end, because R's glm says so out loud. With 25 features and only a few dozen observations the two classes are almost always perfectly separable, and unregularised maximum likelihood has no finite solution there — the coefficients diverge and glm.fit reports that it did not converge. The error rate it produces is whatever IRLS reached before giving up.
That is not an artefact to be hidden: separation is the small-sample failure of a discriminative fit, and it is precisely the variance problem Ng & Jordan describe. But it should be named, and its frequency counted, rather than left in the warning stream.
set.seed(1); p<-25
gen<-function(n){ y<-rbinom(n,1,0.5); X<-matrix(rnorm(n*p),n,p)+y*0.45; data.frame(X, y=factor(y)) }
tst<-gen(3000); ns<-c(30,80,250,800,3000); nbE<-lgE<-numeric(length(ns))
for(k in seq_along(ns)){ a<-b<-c()
for(r in 1:8){ tr2<-gen(ns[k]); if(length(unique(tr2$y))<2) next
pn<-predict(naiveBayes(y~.,tr2), tst); a<-c(a, mean(pn!=tst$y))
g<-glm(y~., tr2, family=binomial); b<-c(b, mean((predict(g,tst,type="response")>=.5)!=(tst$y=="1"))) }
nbE[k]<-mean(a); lgE[k]<-mean(b) }
options(repr.plot.width=8, repr.plot.height=4.6)
plot(ns, nbE, type="b", pch=19, col="#2f855a", lwd=2, log="x", ylim=range(c(nbE,lgE)), xlab="training size n (log)", ylab="test error", main="Ng-Jordan: generative (NB) converges faster")
lines(ns, lgE, type="b", pch=15, col="#2b6cb0", lwd=2)
legend("topright", c("Naive Bayes (generative)","Logistic (discriminative)"), col=c("#2f855a","#2b6cb0"), pch=c(19,15), lwd=2, bty="n")
cat(sprintf("n=%d: NB error %.3f < logistic %.3f; n=%d: logistic %.3f catches NB %.3f.\n", ns[1], nbE[1], lgE[1], ns[length(ns)], lgE[length(ns)], nbE[length(ns)]))
cat("\nHow often was the logistic fit actually a converged one?\n")
for(k in seq_along(ns)){
nsep <- 0
for(r in 1:20){
t3 <- gen(ns[k]); g3 <- suppressWarnings(glm(y~., t3, family=binomial))
if(!g3$converged || max(abs(coef(g3)), na.rm=TRUE) > 50) nsep <- nsep + 1
}
cat(sprintf(" n=%5d: %2d of 20 fits separated or failed to converge\n", ns[k], nsep))
}
cat("\nSeparation is worst in the middle of the small-n range, not at the very bottom: at n=80 almost every fit fails,\n")
cat("while at n=30 fewer do -- with 25 features and 30 rows the design matrix is nearly saturated and the fit often\n")
cat("collapses before separation can be diagnosed. By n=250 the problem is gone entirely.\n")
cat("\nWhere it bites, the logistic error is not a converged estimate: the classes are separable and the likelihood has\n")
cat("no maximum. Naive Bayes has no such failure mode -- it estimates a mean and a variance per feature per class and\n")
cat("those always exist. That IS the Ng-Jordan variance argument in its most extreme form, and it is worth stating\n")
cat("plainly rather than leaving it to a warning message.\n")
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
Warning message: "glm.fit: algorithm did not converge"
Warning message: "glm.fit: fitted probabilities numerically 0 or 1 occurred"
n=30: NB error 0.297 < logistic 0.362; n=3000: logistic 0.126 catches NB 0.126.
How often was the logistic fit actually a converged one?
n= 30: 3 of 20 fits separated or failed to converge n= 80: 19 of 20 fits separated or failed to converge n= 250: 0 of 20 fits separated or failed to converge n= 800: 0 of 20 fits separated or failed to converge n= 3000: 0 of 20 fits separated or failed to converge
Separation is worst in the middle of the small-n range, not at the very bottom: at n=80 almost every fit fails,
while at n=30 fewer do -- with 25 features and 30 rows the design matrix is nearly saturated and the fit often
collapses before separation can be diagnosed. By n=250 the problem is gone entirely.
Where it bites, the logistic error is not a converged estimate: the classes are separable and the likelihood has
no maximum. Naive Bayes has no such failure mode -- it estimates a mean and a variance per feature per class and
those always exist. That IS the Ng-Jordan variance argument in its most extreme form, and it is worth stating
plainly rather than leaving it to a warning message.
3. Reliability and summary¶
The reliability diagram on the credit data confirms the miscalibration: logistic regression tracks the 45° line, Naive Bayes' curve is pushed toward the axes (over-confident 0/1 probabilities). Same good ranking, wrong numbers.
rel<-function(p,y){ y<-as.integer(y)-1; b<-cut(p,seq(0,1,0.1),include.lowest=TRUE)
data.frame(x=tapply(p,b,mean), y=tapply(y,b,mean)) }
rn<-rel(pnb,te$default); rl<-rel(plg,te$default)
options(repr.plot.width=7.5, repr.plot.height=4.6)
plot(c(0,1),c(0,1),type="l",lty=2,xlab="mean predicted probability",ylab="observed default frequency",main="Reliability: NB overconfident, logistic calibrated")
lines(rn$x, rn$y, type="b", pch=19, col="#2f855a", lwd=2); lines(rl$x, rl$y, type="b", pch=15, col="#2b6cb0", lwd=2)
legend("topleft", c(sprintf("Naive Bayes (ECE %.2f)",ece(pnb,te$default)), sprintf("Logistic (ECE %.2f)",ece(plg,te$default))), col=c("#2f855a","#2b6cb0"), pch=c(19,15), lwd=2, bty="n")
Summary¶
e1071, MASS, and glm reproduce the Python notebook's conclusions, and pin down two things the Python side leaves implicit.
First, the two languages' tables differ only because of the split. Given Python's exact train/test indices, e1071 returns its numbers to three decimals. Across twenty random splits Naive Bayes' AUC is stable to a standard deviation of 0.007, while its accuracy swings with a standard deviation of 0.069 — sixteen times logistic's 0.004 — because a model whose probabilities pile up at 0 and 1 has an accuracy that hinges on how many points fall either side of the cut. The ranking metric replicates and the threshold metric does not, which is the miscalibration story in one line.
Second, much of the small-$n$ end of the learning curve rests on logistic fits that are not converged: at $n=80$, nineteen of twenty separate or fail outright, since with 25 features and a few dozen rows the classes are linearly separable and unregularised maximum likelihood has no finite solution. By $n=250$ the problem has vanished. That is not a flaw in the experiment — separation is the discriminative small-sample failure in its most extreme form, and Naive Bayes has no equivalent, since a mean and a variance always exist. It is the Ng-Jordan variance argument taken to its limit, and better said out loud than left in a warning stream. Guidance: generative for small/wide/streaming problems and speed; discriminative when data are plentiful and calibrated probabilities matter (and recalibrate a generative classifier first). This opens the ML arc and bridges the portfolio's Bayesian generative core to its discriminative methods; cross-links to the Calibration and Regularized & Kernel notebooks.