積分大概根據小數位不同

Question

我試圖將牛頓公式應用於長度n = 500的一些反抽樣方法。 （u均勻隨機數，Fphi_x我們根據分佈，fphi_x我們的密度）。 牛頓式則是：求x，使得Fphi_x（X）-u = 0

`u=runif(500) 
xnewton=rep(1e-30,500) 
fphi_x=function(x){return(0.307/(x^{5/12}*(x+1)))} 
Fphi_x=function(u){k=integrate(fphi_x,0,u)$value 
       return(k)} 
for(i in 1:500) 
{ 
    while(abs(Fphi_x(xnewton[i])-u[i])>0.0000001){ 
    xnewton[i]=xnewton[i]-(Fphi_x(xnewton[i])-u[i])/fphi_x(xnewton[i]) 
    } 
}` 

的代碼是正確。但我得到的錯誤代碼： Error in integrate(fphi_x, 0, u) : the integral is probably divergent. 我發現：

Fphi_x(1e15)=1.197390096030634e-05 ，但如果我有它的小數;例如 Fphi_x(10000000000.000000)我收到錯誤消息。 根據我的分佈，我將確保獲得高數字，並帶有小數位。 431232.132131。如何解決這個問題？ 致以問候

Answer 1

這個問題似乎與您的被積函數在下限無限的事實有關;即

> fphi_x(0) 
[1] Inf 

顯然integrate可以處理這個較小的上限但不是較大的上限。我嘗試了幾種方法。一說，似乎工作就是讓變量的變化在整體

 y = 1/x 

其保持積有限，但使得上限明確無窮。然後，您可以使用pracma包中的積分函數quadinf，該函數允許積分上的無限界限。
爲了比較來自兩個集成例程的結果，使用了set.seed(100)，它修復了爲u返回的隨機值。對於u這些值，integrate函數將處理u的前126個值，而quadinf將計算所有500個值。每個的前126個值是相同的。
代碼既如下：

set.seed(100) 
    u=runif(500) 
    xnewton=rep(1e-30,500) 
    fphi_x=function(x){return(0.307/(x^{5/12}*(x+1)))} 
    Fphi_x=function(u){ k=integrate(fphi_x,0,u)$value 
         return(k)} 
    for(i in 1:126) 
    { 
    while(abs(Fphi_x(xnewton[i])-u[i])>0.0000001){ 
     xnewton[i]=xnewton[i]-(Fphi_x(xnewton[i])-u[i])/fphi_x(xnewton[i]) 
    } 
    } 
    xnewton_x <- xnewton 
# 
# using change of variable y = 1/x and quadinf function 
# 
    library(pracma) 
    xnewton=rep(1e-30,500) 
    fphi_y <- function(y) .307*y^(5/12-1)/(1+y) 
    Fphi_y <- function(u) quadinf(fphi_y, 1/u, Inf)$Q 

    for(i in 1:500) 
    { 
    while(abs(Fphi_y(xnewton[i])-u[i])>0.0000001){ 
     xnewton[i] <- xnewton[i]- (Fphi_y(xnewton[i]) - u[i])/fphi_x(xnewton[i]) 
    } 
    } 
    xnewton_y <- xnewton 

quadinf似乎需要超過integrate運行，但至少似乎返回結果情況下，當integrate失敗。