Integral error in R : a limit is NA or NaN-CodePudding

The difference between the two functions below fInt1 and fInt2 is the additional multiplicative term df$ui[i]. Integrating fInt1 works and gives a solution. However, integrating fInt2 gives a limit is NA or NAN error. Where could I be going wrong? This is a follow-up question of Computing integral of a function which is a multiplication of two other integral functions in R

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

Function for the inner integral

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

Function yielding an error

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

CodePudding user response：

The error is the indices in the sapply. It should use 1:nrow(df) instead of seq_along(df), which is 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