So i have this code in R that tries to double integrate a function f(x,y) using integral2, where one domain is (0.5 , +inf), but integral2 doesn't support it (because of that infinite). My question is: do you know a way to make it work?
fprob8 <- function(a , b)
{
f <- function(x , y)
{
densGamma <- function(x)
{
numitor <- b ^ a * fgam(a)
numarator <- x ^ (a - 1) * exp(-x / b)
return (numarator / numitor)
}
densBeta <- function(x)
{
numitor <- fbet(a , b)
numarator <- (1 - x) ^ (b - 1) * x ^ (a - 1)
return (numarator/numitor)
}
return (densGamma(x) * densBeta(y))
}
ymin <- function(x) # limita superioara pentru cea de-a doua integrala
{
return (min(x - 0.5),1)
}
I <- integral2(f, 0.5 , Inf , ymin , 1)
return(I)
}
Use the cubature
package. It evaluates multiple integrals, allowing infinite bounds.
Here is an example:
library(cubature)
f <- function(x) exp(-x[1]-x[2])
pcubature(f, c(0,0), c(Inf,Inf))$integral
# 1
library(pracma)
f2 <- function(x,y) f(c(x,y))
integral2(f2, 0, 1000, 0, 1000, vectorized = FALSE)$Q
# 1