嵌套Logit模型序贯估计中的奇异性问题排查求助

错误原因分析 #

你遇到的奇异性错误（Lapack routine dgesv: system is exactly singular）核心原因是参数不可识别：

• 你的数据中，每个观测只属于一个巢：要么是beach/pier（岸线巢），要么是boat/charter（船只巢），没有观测同时包含两个巢的选项。
• 当你运行mlogit(mode ~ price, fish.joint)时，模型默认会尝试估计所有四个选项的截距（以其中一个为基准），但由于两个巢的观测完全互斥，岸线巢的截距参数和船只巢的截距参数之间没有重叠的观测数据支撑，导致设计矩阵出现线性相关性，Hessian矩阵无法求逆，最终触发奇异性错误。

解决办法 #

方法1：约束基准选项截距，实现可识别的联合估计 #

为每个巢设定一个基准选项（比如岸线巢用beach，船只巢用boat），将它们的截距固定为0，只估计巢内其他选项的相对截距，同时共享价格系数。这样参数就变得可识别了：

library(mlogit)
data("Fishing", package = "mlogit")

# 你的原始数据处理步骤
Fishing$avail.beach <- 0
Fishing$avail.pier <- 0
Fishing$avail.boat <- 0
Fishing$avail.charter <- 0
Fishing$avail.charter[Fishing$mode == "charter" | Fishing$mode == "boat"] <- 1
Fishing$avail.boat[Fishing$mode == "charter" | Fishing$mode == "boat"] <- 1
Fishing$avail.pier[Fishing$mode == "pier" | Fishing$mode == "beach"] <- 1
Fishing$avail.beach[Fishing$mode == "pier" | Fishing$mode == "beach"] <- 1

Fish <- mlogit.data(Fishing, shape="wide", varying=c(2:9, 11:14), choice="mode")
fish.joint <- subset(Fish, avail ==1)

# 约束基准选项的截距为0，解决参数不可识别问题
const_pars <- c("beach:(intercept)" = 0, "boat:(intercept)" = 0)
fish.1 <- mlogit(mode ~ price, data = fish.joint, constPar = const_pars)

# 查看估计结果
summary(fish.1)

方法2：直接使用嵌套Logit模型原生设定 #

如果你最终目标是估计嵌套Logit模型，无需手动拆分可用选项，直接用mlogit的nests参数指定嵌套结构，函数会自动处理巢内和巢间的参数估计，避免奇异性：

library(mlogit)
data("Fishing", package = "mlogit")

# 整理宽格式数据
Fish <- mlogit.data(Fishing, shape="wide", varying=c(2:9), choice="mode")

# 定义嵌套结构
nests_list <- list(shore = c("beach", "pier"), boat = c("boat", "charter"))

# 估计嵌套Logit模型
fish_nl <- mlogit(mode ~ price, data = Fish, nests = nests_list, un.nest.el = TRUE)

# 查看结果
summary(fish_nl)

方法3：分别估计子模型再联合（进阶） #

如果你坚持按Ben Akiva & Lerman的步骤手动联合估计，可以先分别估计每个巢的MNL模型，提取参数作为联合估计的初始值，再通过参数约束实现联合估计：

# 分别估计两个子模型
model_shore <- mlogit(mode ~ price, data = subset(fish.joint, alt %in% c("beach", "pier")))
model_boat <- mlogit(mode ~ price, data = subset(fish.joint, alt %in% c("boat", "charter")))

# 提取初始参数
start_vals <- c(coef(model_shore), coef(model_boat))

# 联合估计（共享价格系数）
fish_joint_manual <- mlogit(mode ~ price | alt, data = fish.joint, 
                            constPar = c("beach:(intercept)"=0, "boat:(intercept)"=0),
                            start = start_vals)

内容的提问来源于stack exchange，提问作者Sven Müller