为何MASS::lda()的类别后验概率与自行计算结果不一致？

问题背景 #

希望在线性判别分析（LDA）框架下计算后验概率，根据给定X选择后验概率最大的类别。使用MASS::lda()得到结果后，自行编写代码计算无法得到一致结果，可复现代码如下：

# x is a numeric vector (variable)  
x <- matrix(c(5,6,5,7,2,1,3,2,4,2,4 ), ncol=1, byrow=T)

# y are the k (classes)
y <- c(1,1,1,1,1,1,1,2,2,2,2)

# merge into data.frame
df <- data.frame(x,y)

# plot values over a line
plot(x, col =  y, pch= y, cex = 1.5)

# run lda with MASS
library(MASS)
mlda <- lda(y ~ ., data = df)

# predict class
m_prev <- predict(mlda)

# see posterior probabilities
m_prev$posterior

# now i run code to compute 
# the posteriors:

# prior mean each k
prior1 <- sum(df$y==1)/nrow(df)
prior2 <- sum(df$y==2)/nrow(df)

# mean(x) each k
m1 <- round(mean(df[df$y==1,"x"]),2)
m2 <- round(mean(df[df$y==2,"x"]),2)

# standard deviation
s = sd(df$x)

# f(x) each k 
fnorm1 <- dnorm(x, mean = m1, sd = s)
fnorm2 <- dnorm(x, mean = m2, sd = s)

# and finaly the posterior probabilites k=1 using the bayes theorem
# p(y=k|x=x) = fk(x)p(k) / sum(fi(x)p(i))
post1 <- prior1 * fnorm1 / (prior1 * fnorm1 + prior2 * fnorm2)

# it is not the same as the posterior probabilities 
# obtained with MASS::lda()
cbind(post1, m_prev$posterior[,1])

错误点分析 #

1. 类内标准差计算错误：LDA假设各类别的类内方差相同，应使用合并类内方差（pooled within-class variance），而非全局标准差。全局标准差包含了类间差异，不符合LDA的假设前提。
2. 均值精度损失：对类均值做了round()操作，引入了不必要的精度误差，导致正态密度函数计算结果偏差。
3. 方差计算细节：合并类内方差的计算需要基于每个类的样本方差，按类内样本量加权平均后再开平方得到标准差。

修正后的代码 #

# 重新计算类均值（不做round）
m1 <- mean(df[df$y==1,"x"])
m2 <- mean(df[df$y==2,"x"])

# 计算合并类内方差
n1 <- sum(df$y==1)
n2 <- sum(df$y==2)
var1 <- var(df[df$y==1,"x"])
var2 <- var(df[df$y==2,"x"])
pooled_var <- ((n1-1)*var1 + (n2-1)*var2)/(n1+n2-2)
s <- sqrt(pooled_var)

# 重新计算正态密度
fnorm1 <- dnorm(x, mean = m1, sd = s)
fnorm2 <- dnorm(x, mean = m2, sd = s)

# 计算后验概率
post1_corrected <- prior1 * fnorm1 / (prior1 * fnorm1 + prior2 * fnorm2)

# 对比结果
cbind(post1_corrected, m_prev$posterior[,1])

结果验证 #

运行修正后的代码后，post1_corrected与m_prev$posterior[,1]的结果会完全一致，解决了之前的差异问题。

内容的提问来源于stack exchange，提问作者maluicr