R 中的积分误差:极限是 NA 或 NaN

Integral error in R : a limit is NA or NaN

fInt1fInt2下面两个函数的区别在于增加了乘法项df$ui[i]。集成 fInt1 有效并给出了解决方案。但是,积分 fInt2 得到 limit is NA or NAN error。我哪里会出错?这是

的后续问题
set.seed(1234)
G <- 5# Suppose 5 groups
theta<-0.5
n_i <- 2 # There are two individuals per group
nTot <- n_i*G # In total we have 4 individuals 
z_ij <- rnorm(nTot, 0, 0.1)
ui<-rnorm(nTot, 0, 0.2)
T_ij <- runif(nTot, 0, 15)
Data <- round(data.frame(id = rep(1:nTot), group = rep(1:G, rep(2,G)), ui,z_ij, T_ij=round(T_ij,1)) , 3)
head(Data)
  id group     ui   z_ij T_ij
1  1     1 -0.095 -0.121  8.3
2  2     1 -0.200  0.028  9.7
3  3     2 -0.155  0.108  4.7
4  4     2  0.013 -0.235  9.3
5  5     3  0.192  0.043  4.9
6  6     3 -0.022  0.051  7.5

内积分函数

fInt1 <- function(df) {
  Vectorize({function(y) {
    prod(
      sapply(
        seq_along(df),
        function(i) integrate(function(x) x*y*exp(x + y + df$z_ij[i]), 0, df$T_ij[i])$value
      )
    )
  }})
}

fInt2 <- function(df) {
  Vectorize({function(y) {
    prod(
      sapply(
        seq_along(df),
        function(i) integrate(function(x) x*df$ui[i]*y*exp(x + y + df$z_ij[i]), 0, df$T_ij[i])$value
      )
    )
  }})
}
GroupInt1 <- sapply(1:G, function(grp) integrate(fInt1(subset(Data, group == grp, select = c("z_ij","T_ij"))), -5, 5)$value)
GroupInt1
[1] [1] 8.579064e+14 7.361849e+12 1.529633e+12 4.659699e+14 2.230921e+13

函数产生错误

GroupInt2 <- sapply(1:G, function(grp) integrate(fInt2(subset(Data, group == grp, select = c("z_ij","ui", "T_ij"))), -5, 5)$value)
Error in integrate(function(x) x * df$ui[i] * y * exp(x + y + df$z_ij[i]), : a limit is NA or NaN

错误是 sapply 中的索引。它应该使用 1:nrow(df) 而不是 seq_along(df),即 column-wise.

set.seed(1234)
G <- 5# Suppose 5 groups
n_i <- 2 # There are two individuals per group
nTot <- n_i*G # In total we have 4 individuals 
z_ij <- rnorm(nTot, 0, 0.1)
ui <- rnorm(nTot, 0, 0.2)
T_ij <- runif(nTot, 0, 15)
Data <- round(data.frame(id = rep(1:nTot), group = rep(1:G, rep(2,G)), ui, z_ij, T_ij=round(T_ij,1)), 3)

fInt1 <- function(df) {
  Vectorize({function(y) {
    prod(
      sapply(
        1:nrow(df),
        function(i) y*exp(y + df$z_ij[i])*integrate(function(x) x*exp(x), 0, df$T_ij[i])$value
      )
    )
  }})
}

fInt2 <- function(df) {
  Vectorize({function(y) {
    prod(
      sapply(
        1:nrow(df),
        function(i) df$ui[i]*y*exp(y + df$z_ij[i])*integrate(function(x) x*exp(x), 0, df$T_ij[i])$value
      )
    )
  }})
}

GroupInt1 <- sapply(1:G, function(grp) integrate(fInt1(subset(Data, group == grp, select = c("z_ij","T_ij"))), -5, 5)$value)
GroupInt1
#> [1] 8.579064e+14 7.361849e+12 1.529633e+12 4.659699e+14 2.230921e+13

GroupInt2 <- sapply(1:G, function(grp) integrate(fInt2(subset(Data, group == grp, select = c("z_ij", "ui", "T_ij"))), -5, 5)$value)
GroupInt2
#> [1]  1.630022e+13 -1.483413e+10 -6.461171e+09  8.650265e+12 -1.799484e+12