数值计算二阶导数误差大
Large error when numerically computing second derivative
我正在尝试在 R 中对函数 l
(泊松分布在向量 x
和 lambda=6
上的对数)进行数值计算,这是我的代码:
x=c(2,3)
t=6
delta=1e-12
h=1e-12
L=function(x,t) dpois(x,t)
l<-function(x,t) log(prod(L(x,t)))
ld<-function(x,t) (l(x,t+delta)-l(x,t))/delta
ldd<-function(x,t) (ld(x,t+h)-ld(x,t))/h
ld(x,t)
ldd(x,t)
我的输出是
> ld(x,t)
[1] -1.167066
> ldd(x,t)
[1] 888178420
但是对于这个完全相同的函数,我使用 wolfram 并得到 -7/6~~-1.16667 for the first derivative and -5/36~~-0.138889 作为二阶导数。在过去的两个小时里,我一直在试图弄清楚为什么我的函数有这么大的错误。
注意:这是一个 class 项目,所以我不能在 R 中使用导数函数。
我认为您看到的问题是由于四舍五入和其他数值问题造成的。我建议两种方法:
- 提高您的
delta
和 h
值;由于您的计算机精度有限,因此在许多情况下 t
和 t+1e-12
将完全相同。当然,这在数值计算导数时会造成很大的问题。
- 为了数值稳定性,将
log(prod(x))
替换为 sum(log(x))
通常是个好主意。
通过这两个调整,我们得到了更好的结果:
delta = 1e-5
h = 1e-4
L=function(x,t) dpois(x,t)
l<-function(x,t) sum(log(L(x,t)))
ld<-function(x,t) (l(x,t+delta)-l(x,t))/delta
ldd<-function(x,t) (ld(x,t+h)-ld(x,t))/h
ld(x,t)
# [1] -1.166667
ldd(x,t)
# [1] -0.1388853
我正在尝试在 R 中对函数 l
(泊松分布在向量 x
和 lambda=6
上的对数)进行数值计算,这是我的代码:
x=c(2,3)
t=6
delta=1e-12
h=1e-12
L=function(x,t) dpois(x,t)
l<-function(x,t) log(prod(L(x,t)))
ld<-function(x,t) (l(x,t+delta)-l(x,t))/delta
ldd<-function(x,t) (ld(x,t+h)-ld(x,t))/h
ld(x,t)
ldd(x,t)
我的输出是
> ld(x,t)
[1] -1.167066
> ldd(x,t)
[1] 888178420
但是对于这个完全相同的函数,我使用 wolfram 并得到 -7/6~~-1.16667 for the first derivative and -5/36~~-0.138889 作为二阶导数。在过去的两个小时里,我一直在试图弄清楚为什么我的函数有这么大的错误。
注意:这是一个 class 项目,所以我不能在 R 中使用导数函数。
我认为您看到的问题是由于四舍五入和其他数值问题造成的。我建议两种方法:
- 提高您的
delta
和h
值;由于您的计算机精度有限,因此在许多情况下t
和t+1e-12
将完全相同。当然,这在数值计算导数时会造成很大的问题。 - 为了数值稳定性,将
log(prod(x))
替换为sum(log(x))
通常是个好主意。
通过这两个调整,我们得到了更好的结果:
delta = 1e-5
h = 1e-4
L=function(x,t) dpois(x,t)
l<-function(x,t) sum(log(L(x,t)))
ld<-function(x,t) (l(x,t+delta)-l(x,t))/delta
ldd<-function(x,t) (ld(x,t+h)-ld(x,t))/h
ld(x,t)
# [1] -1.166667
ldd(x,t)
# [1] -0.1388853