使用R语言nls拟合非线性指数模型时遇错误的技术求助

问题说明 #

需要为变量y和x构建形式为(y = a \cdot e^{r \cdot x})的非线性指数模型并提取残差，但使用R的nls函数时遇到两类错误：

• 直接建模报错：singular gradient（奇异梯度）
• 指定alg="plinear"后报错：step factor 0.000488281 reduced below 'minFactor' of 0.000976562

数据集如下：

df <- structure(list(y = c(464208.56, 334962.43, 361295.68, 426535.68, 258843.93, 272855.46, 
   166322.72, 244695.28, 227003.03, 190728.4, 156025.45, 72594.24, 56911.4, 175328.95, 161199.76, 
   152520.77, 190610.57, 60734.34, 31620.9, 74518.86, 45524.49, 2950.58, 2986.38, 15961.77, 12484.05, 
   6828.41, 2511.72, 1656.12, 5271.4, 7550.66, 3357.71, 3620.43, 3699.85, 3337.56, 4106.55, 3526.66, 
   2996.79, 1649.89, 4561.64, 1724.25, 3877.2, 4426.69, 8557.61, 6021.61, 6074.17, 4072.77, 4032.95, 
   5280.16, 7127.22), 
   x = c(39.23, 38.89, 38.63, 38.44, 38.32, 38.27, 38.3, 38.4, 38.56, 38.79, 39.06, 39.36, 39.68, 
   40.01, 40.34, 40.68, 41.05, 41.46, 41.93, 42.48, 43.14, 43.92, 44.84, 45.9, 47.1, 48.4, 49.78, 
   51.2, 52.62, 54.01, 55.31, 56.52, 57.6, 58.54, 59.33, 59.98, 60.46, 60.78, 60.94, 60.92, 60.71, 
   60.3, 59.69, 58.87, 57.86, 56.67, 55.33, 53.87, 52.33)), 
   row.names = c(NA, -49L), 
   class = c("tbl_df", "tbl", "data.frame"), 
   na.action = structure(c(`1` = 1L, `51` = 51L), 
   class = "omit"))

解决方法 #

1. 优化初始值（核心解决思路） #

nls对初始参数值极度敏感，原初始值a=0.5与实际数据量级差距过大（y最大值达46万，按原初始值计算的预测值远小于实际值），导致算法无法收敛。可以通过线性化转换先估算合理初始值：

# 对指数模型取对数，转换为线性模型估计初始参数
lm_fit <- lm(log(y) ~ x, data = df)
init_a <- exp(coef(lm_fit)[[1]])  # 将线性模型的截距转换为a的初始值
init_r <- coef(lm_fit)[[2]]       # 线性模型的斜率即为r的初始值

# 使用优化后的初始值拟合nls模型
m <- nls(y ~ a * exp(r * x), 
         start = list(a = init_a, r = init_r), 
         data = df)

# 查看模型拟合结果
summary(m)

# 提取残差
residuals(m)

2. 稳健拟合备选方案 #

如果优化初始值后仍有问题，可使用minpack.lm包的nlsLM函数，它采用Levenberg-Marquardt算法，对初始值要求更低，拟合稳定性更强：

# 安装并加载包
install.packages("minpack.lm")
library(minpack.lm)

# 拟合模型
m_lm <- nlsLM(y ~ a * exp(r * x), 
              start = list(a = init_a, r = init_r), 
              data = df)

# 查看结果及残差
summary(m_lm)
residuals(m_lm)

3. 模型效果验证 #

拟合后可通过绘图直观验证模型拟合效果：

plot(y ~ x, data = df, main = "指数模型拟合结果", pch = 16)
lines(df$x, predict(m, newdata = df), col = "red", lwd = 2)

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