1

通常の逆ガウス最尤推定を行います。optim()私が得た後

「エラー×統合（f2、下限= 0、上限= Inf、z = z）：関数の評価で間違った長さの結果が得られました」.

バグのあるコードの部分は

f2 <- function (t, z) {
    exp(-t)*t^(1-1/2)*(1+t/2*z)^(1-1/2)
}
int2 <- function (z) { 
    integrate(f2,  lower=0, upper=Inf, z=z)$value 
}

R サイトでの同様の推奨事項には、まだ同じエラーがあります。

int2 <- function (z) { 
    integrate(f2,  lower=0, upper=Inf, z=z) 
}

問題は、R のパラメーター z (数値ではなく記号) に応じて、t に関して f2 を統合する方法です。

コード全体は以下

#NIG
library(stats) 
unibubble_nig_ur <- function(x){
    #mu
    x1 <- x[1]
    #nu
    x2 <- x[2]
    #alpha
    x3 <- x[3]
    #beta
    x4 <- x[4]
    #kappa
    k <- x[5]
    #sigma_sq
    x6 <- x[6]
    H <- log((x3^x4+t2^x4)/(x3^x4+t1^x4))
    #mu_t
    x7 <- x1+(x2-x2^2/2)*H*k
    #alpha_t
    x8 <- sqrt(1/((x6-x2^2*H))*k+1/4)
    #beta_t
    x9<--c(1/2)
    #delta_t
    x10 <- sqrt((x6-x2^2*H)/k)
    zval <- x8*sqrt(x10^2+(yobs-x7^2))

    f = function (t) exp(-t)*(1)^(t-1)
    gamma = integrate(f, lower=0, upper=Inf)

    #Bessel function of 2nd kind
    f2 <- function (t, z) {
        exp(-t)*t^(1-1/2)*(1+t/2*z)^(1-1/2)
    }
    int2 <- function (z) { 
        integrate(f2,  lower=0, upper=Inf, z=z)$value 
    }

    K1 <- exp(-zval)*sqrt(pi/2*zval)*int2(zval)/gamma
    n <- length(logprice)
    t1 <- seq(1, n-1)
    t2 <- 1+t1
    yobs <- diff(logprice,1)

    # Probability Density Function 
    dnig <- x8*x10*exp(x10*sqrt(x8^2-x9^2)+x9*(yobs-x7))/(pi*sqrt(x10^2+(yobs-x7^2)))
    loglklhd<--sum(log(dnig))
    loglklhd
}

out_nig_ur <- optim(unibubble_nig_ur, p=c(60,40,80,80,80,80), hessian=TRUE)
out_nig_ur
4

3 に答える 3

4

例えば：

F2 <- function(z){
    f2 <- function(t){
       exp(-t)*t^(1-1/2)*(1+t/2*z)^(1-1/2)
   }
return(f2)
}
int2 <- function (z) { 
    integrate(F2(z),  lower=0, upper=Inf)$value 
}
> int2(1)
[1] 1.156935
于 2012-06-24T08:35:31.403 に答える
1

ベクトル引数にint2を適用するには、新しい関数を導入して、数式int2vecの代わりにそれを使用できます。int2

int2vec<-function(t) sapply(t,int2)
于 2012-06-22T18:49:20.557 に答える
1

yobs長さが 1 より大きい場合も同様であり、行内のへの呼び出しzvalの戻り値も同様であり、混乱します。f2int2K1 <- exp(-zval)*sqrt(pi/2*zval)*int2(zval)/gammaintegrate

于 2012-06-22T14:17:34.703 に答える