r - 有没有比 uniroot 更精确的求解器来求解非线性方程？

Question

我的问题是我使用 uniroot 来求解非线性函数的值，但问题是当我手动计算时这些值与预期不符。

我写了一个函数，如：

yield.to.call <- function(ytc,p,c,t,f=4,N=100){
  #d <- (1+(ytc/f))^f
  c <- c/f 
  if(t<92){t <- ((92-t)/(92-0))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+c/((1+(ytc/f))^(t+13))+c/((1+(ytc/f))^(t+14))+c/((1+(ytc/f))^(t+15))+c/((1+(ytc/f))^(t+16))+c/((1+(ytc/f))^(t+17))+c/((1+(ytc/f))^(t+18))+(c+N)/((1+(ytc/f))^(t+19))-p}}
  else if(t<184){t <- ((184-t)/(184-92))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+c/((1+(ytc/f))^(t+13))+c/((1+(ytc/f))^(t+14))+c/((1+(ytc/f))^(t+15))+c/((1+(ytc/f))^(t+16))+c/((1+(ytc/f))^(t+17))+(c+N)/((1+(ytc/f))^(t+18))-p}}
  else if(t<275){t <- ((275-t)/(275-184))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+c/((1+(ytc/f))^(t+13))+c/((1+(ytc/f))^(t+14))+c/((1+(ytc/f))^(t+15))+c/((1+(ytc/f))^(t+16))+(c+N)/((1+(ytc/f))^(t+17))-p}}
  else if(t<365){t <- ((365-t)/(365-275))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+c/((1+(ytc/f))^(t+13))+c/((1+(ytc/f))^(t+14))+c/((1+(ytc/f))^(t+15))+(c+N)/((1+(ytc/f))^(t+16))-p}}
  else if(t<457){t <- ((457-t)/(457-365))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+c/((1+(ytc/f))^(t+13))+c/((1+(ytc/f))^(t+14))+(c+N)/((1+(ytc/f))^(t+15))-p}}
  else if(t<548){t <- ((548-t)/(548-457))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+c/((1+(ytc/f))^(t+13))+(c+N)/((1+(ytc/f))^(t+14))-p}}
  else if(t<639){t <- ((639-t)/(639-548))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+c/((1+(ytc/f))^(t+12))+(c+N)/((1+(ytc/f))^(t+13))-p}}
  else if(t<730){t <- ((730-t)/(730-639))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+c/((1+(ytc/f))^(t+11))+(c+N)/((1+(ytc/f))^(t+12))-p}}
  else if(t<821){t <- ((821-t)/(821-730))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+c/((1+(ytc/f))^(t+10))+(c+N)/((1+(ytc/f))^(t+11))-p}}
  else if(t<915){t <- ((915-t)/(915-821))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+c/((1+(ytc/f))^(t+9))+(c+N)/((1+(ytc/f))^(t+10))-p}}
  else if(t<1006){t <- ((1006-t)/(1006-915))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+c/((1+(ytc/f))^(t+8))+(c+N)/((1+(ytc/f))^(t+9))-p}}
  else if(t<1097){t <- ((1097-t)/(1097-1006))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+c/((1+(ytc/f))^(t+7))+(c+N)/((1+(ytc/f))^(t+8))-p}}
  else if(t<1188){t <- ((1188-t)/(1188-1097))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+c/((1+(ytc/f))^(t+6))+(c+N)/((1+(ytc/f))^(t+7))-p}}
  else if(t<1280){t <- ((1280-t)/(1280-1188))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+c/((1+(ytc/f))^(t+5))+(c+N)/((1+(ytc/f))^(t+6))-p}}
  else if(t<1371){t <- ((1371-t)/(1371-1280))#done 
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+c/((1+(ytc/f))^(t+4))+(c+N)/((1+(ytc/f))^(t+5))-p}}
  else if(t<1461){t <- ((1461-t)/(1461-1371))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+c/((1+(ytc/f))^(t+3))+(c+N)/((1+(ytc/f))^(t+4))-p}}
  else if(t<1553){t <- ((1553-t)/(1553-1461)) 
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+c/((1+(ytc/f))^(t+2))+(c+N)/((1+(ytc/f))^(t+3))-p}}
  else if(t<1645){t <- ((1645-t)/(1645-1553))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+c/((1+(ytc/f))^(t+1))+(c+N)/((1+(ytc/f))^(t+2))-p}}
  else if(t<1736){t <- ((1736-t)/(1736-1645))
    q <- function(ytc,p,c,t,f=4,N=100){c/((1+(ytc/f))^(t))+(N+c)/((1+(ytc/f))^(t+1))-p}}
  else if(t<1826){t<-((1826-t)/(1826-1736))
    q <- function(ytc,p,c,t,f=4,N=100){(N+c)/(1+(ytc/f)^(t))-p}}
  else q <- function(ytc,p,c,t=0,f=4,N=100){(N+c)/p-(ytc/f)-1} 
  
  q2<- function(ytc){q(ytc,p,c,t,N=100)}
  yield.t.c <- uniroot(q2,c(0,1),tol = .Machine$double.eps^0.5)$root 
  return(yield.t.c)
}

yield.to.call(ytc,p=100,t=1790,c=5,N=100)

例如 t=1790 的值是 6.987728e-05 但应该是 0.0126 左右。

调用收益率是一个函数，它应该为我计算现金流 a 息票债券的内部收益率。意思是结构 y=c/(1+(r/4))^t 的现金流量，其中除了我试图解决的 r 之外的一切都是已知的。

有人有想法吗？

编辑：我发现问题出在以下行：

else if(t<1826){t<-((1826-t)/(1826-1736))
    q <- function(ytc,p,c,t,f=4,N=100){(N+c)/(1+(ytc/f)^(t))-p}}

由于函数对于任何 t 值不再是线性的：

t=1790 的图

知道如何“线性化”它并获得我想要的解决方案吗？

标签： rsolver