非线性回归嵌套迭代程序报错：需数值/复数矩阵/向量参数

模型与迭代方案 #

从实证数据获取时间序列（chronosequences），为同步计算**公式(1)**的立地参数与全局参数，采用嵌套迭代程序：

公式(1)：
$$H_2 = H_1 \times \left[ \frac{1 - \exp(b_2 \cdot T_2)}{1 - \exp(b_2 \cdot T_1)} \right]^{(b_3 \cdot (H_1)^{b_1})}$$

参数说明：

• $H_2$：$T_2$龄的实测树高
• $H_1$：$T_1$龄林分高度的立地参数
• $b_1、b_2、b_3$：待估全局参数

嵌套迭代步骤：

1. 以实证数据确定的$H_1$初始值校准全局参数
2. 将全局参数初始值视为常数，估算每株树木的立地参数
3. 以$H_1$估算值为常数重新拟合全局参数，重复迭代至参数稳定

参数估计采用R语言的nls（非线性最小二乘法）程序。

实现代码 #

###### data包含从实证时间序列数据中提取的树木树高与年龄信息
fit1 <- nls(
  H2 ~ H1 * ((1 - exp(b2*T2))/(1 - exp(b2*100)))^(b3*(H1^b1)), 
  data = data,    #### 对应非线性公式(1)
  start = list(b1=-1, b2=-0.04, b3=53),
  control = nls.control(maxiter = 1000, tol=1e-7)
)   
fit1  #### 查看首次拟合参数

# 首次拟合得到的参数值
b1 = -0.8709, b2 = -0.02852, b3 = 28.61296  

predH <- predict(fit1, level = c(0:2))   ######## 第一轮迭代的拟合值
modDat <- data.frame(data, predH)       ######## 将拟合值合并到原数据
modDat ##### 生成新数据框

para <- c(-0.8709, -0.02852, 28.61296)    
para
iter <- 0
repeat {
  b1 <- para[1]
  b2 <- para[2]
  b3 <- para[3]

  Db <- cbind(
    modDat$H1 * ((1 - exp(b2*modDat$T2))/(1 - exp(b2*100)))^(b3*(modDat$H1^b1)),
    modDat$predH * ((1 - exp(b2*modDat$T2))/(1 - exp(b2*100)))^(b3*(modDat$predH^b1))
  )
  new.fit <- para + Solve((Db) %*% as.numeric(sapply(Db, as.numeric)) %*% (Db)) %*% 
    (predH * ((1 - exp(b2$T2))/(1 - exp(b2*100)))^(b3*(predH^b1)))
  iter <- iter + 1
  if (sum((para - new.fit)^2) < 1e-5 || iter > 1000) break
  para <- new.fit
}

# 需先加载systemfit包以使用Solve函数
library(systemfit)

报错问题 #

已尝试用as.numeric转换参数，但仍触发报错：

"requires numeric/complex matrix/vector arguments"

请求解决该问题。

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