使用 nls - 奇异梯度将 erf 拟合到 R 中的几个实验数据点

Fitting erf to few experimental data points in R with nls - singular gradient

提问人:ikempf 提问时间:11/14/2023 最后编辑:ikempf 更新时间:11/15/2023 访问量:35

问:

我需要拟合一个误差函数,该函数仅将数据的物理模型描述为 6 个实验数据点。

error 函数为:

func_erf <- function(x, #m
                 D,     #m2/s
                 t,     #s
                 s      #m
){ 
  result = (erf((x+s)/sqrt(8*D*t)) - erf((x-s)/sqrt(8*D*t)))/(2*erf(s/sqrt(8*D*t)))
  return(result)
}

我试图拟合的数据是:

> data_exp
1 1.000000000 0.000
2 0.766766619 0.001
3 0.252337795 0.002
4 0.098405369 0.003
5 0.046523446 0.004
6 0.004363998 0.005
> dput(data_exp)
structure(list(y = c(1.00000000046026, 0.766766619156469, 0.252337794969704, 
0.0984053685324868, 0.0465234458835242, 0.00436399807604814), 
    x = c(0, 0.001, 0.002, 0.003, 0.004, 0.005)), class = "data.frame", row.names = c(NA, 
-6L))

通过实验,我们知道 t = 6.504601e-05 秒。因此,我们只需要将参数 D 和 s 拟合到我们的实验数据中。

假设起始参数与曲线拟合有点接近实验数据,并用达塔点绘制初始拟合猜测,我得到:enter image description here

然而,nls 拟合过程总是会导致奇异梯度矩阵的误差。

coef_fit_erf_guess = c(1e-8,       #D, #m2/s
                       0.75*1e-3   #s, #m
)
t_exp = 6.504601e-05 #seconds

fit_nls_erf<- nls(y~func_erf(x,D,t= t_exp, s),data=data_exp,
                             start=list(D = coef_fit_erf_guess [1],
                                        s = coef_fit_erf_guess [2]))


为什么这个拟合程序不起作用?我能改进什么?有没有办法找到一个更好的拟合猜测或以迭代的“手动”方式拟合它?

非常感谢您的帮助!

R 曲线拟合 NLS 收敛

评论

1赞 Roland 11/14/2023
我成功地拟合了 ,但我不相信它找到了全局最小值。这些系数似乎具有很强的相关性。minpack.lm::nlsLMcoef_fit_erf_guess = c(0.006, 1e-6)
1赞 user2554330 11/14/2023
您应该使用,以便其他人更容易输入您的数据。dput(data_exp)
0赞 ikempf 11/15/2023
谢谢@user2554330,我已经相应地更新了问题。

答:

0赞 user2554330 11/14/2023 #1

错误表面看起来非常奇怪。您可以使用以下代码获得可旋转的图:

data_exp <- structure(list(id = c(1, 2, 3, 4, 5, 6), 
  y = c(1, 0.766766619,0.252337795, 0.098405369, 0.046523446, 0.004363998), 
  x = c(0, 0.001, 0.002, 0.003, 0.004, 0.005)), 
  class = "data.frame", row.names = c(NA, -6L))

func_erf <- function(x, #m
                     D,     #m2/s
                     t,     #s
                     s      #m
){ 
  result = (erf((x+s)/sqrt(8*D*t)) - erf((x-s)/sqrt(8*D*t)))/(2*erf(s/sqrt(8*D*t)))
  return(result)
}

SS <- function(D, s) {
  with(data_exp, sum((y - func_erf(x, D, t = 6.504601e-05, s))^2))
}

SS <- Vectorize(SS)
library(rgl)
library(pracma)
persp3d(SS, xlim = c(0, 1e-5), ylim = c(0.00095, 0.003))

创建于 2023-11-14 with reprex v2.0.2

我对此感到不惊讶。nls()

评论

0赞 ikempf 11/15/2023
感谢您的回复!您期望哪种类型的误差面才能实现良好的拟合?我同意它充其量看起来很“特别”:)
1赞 user2554330 11/15/2023
最简单的最大化是抛物线曲面。你可以通过线性模型来获得它。一些非线性模型给出的曲面近似为抛物线,而您的曲面可能位于最优值的足够近的邻域,但全局形状与抛物线相去甚远。