Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- # Marja de eroare pentru comparatii pe numere reale.
- EPSILON <- 10e-9
- CONST_DISP <- 5
- # De cate ori verificam sa fie PDF pe fiecare interval al functiei date.
- TESTS <- 500
- equalToOne <- function (x)
- {
- return (x < 1 + EPSILON & x > 1 - EPSILON)
- }
- calcIntegralSum <- function (myFunc, domain)
- {
- sum <- 0
- for (i in domain)
- sum <- sum + integrate (Vectorize (myFunc), i[1], i[2])$value
- return (sum)
- }
- # Ca o functie sa fie densitate de probabilitate trebuie sa fie pozitiva, iar integrala ei pe R sa fie 1.
- checkDensProb <- function (myFunc, domain = list (c (-Inf, Inf)), constNorm = 1)
- {
- tryCatch (
- {
- # Calculam suma integralelor, similar cu cerinta 1.
- integralSum <- calcIntegralSum (myFunc, domain)
- # Comparam cu 1 folosindu-ne de EPSILON.
- if (equalToOne (integralSum / constNorm))
- {
- # Pentru fiecare interval al functiei, verificam daca este PDF.
- for (i in domain)
- {
- # Verificam intervalul [st, dr].
- st <- max(i[1], -2e9)
- dr <- min(i[2], 2e9)
- length <- dr / CONST_DISP - st / CONST_DISP
- nrs <- runif (TESTS, st, dr-length)
- # Testam TESTS subintervale de lungime egala.
- for (mn in nrs)
- {
- mx = mn + length
- # Daca integrala este negativa, functia nu este PDF.
- if (integrate (Vectorize (myFunc), mn, mx)$value < 0) return (FALSE)
- }
- }
- return (TRUE)
- }
- return (FALSE)
- },
- error = function (e)
- {
- # Daca integrarea a esuat, returnam false.
- return (FALSE)
- })
- }
- isPdf <- function(f, suport = list(c(-Inf, Inf)), normCt = 1, TCNT = 1000){
- tryCatch({
- integratedVal <- 0
- for(interval in suport){
- integratedVal <- integratedVal + integrate(Vectorize(f), lower = interval[1], upper = interval[2])$value
- }
- if((integratedVal / normCt) < 1 + 10e-9 & (integratedVal / normCt) > 1 - 10e-9){
- cnt <- TCNT / length(suport)
- for(interval in suport){
- lowRand <- 0
- upperRand <- 0
- if(interval[1] == -Inf){
- lowRand <- -.Machine$integer.max
- }
- else{
- lowRand <- interval[1]
- }
- if(interval[2] == Inf){
- upperRand <- .Machine$integer.max
- }
- else{
- upperRand <- interval[2]
- }
- testLen <- upperRand / 10 - lowRand / 10
- rands <- runif(cnt, lowRand, upperRand - testLen)
- for(x in c(1:cnt)){
- lowerBound <- rands[x]
- upperBound = lowerBound + testLen
- if(lowerBound > upperBound){
- swapAux <- upperBound
- upperBound <- lowerBound
- lowerBoud <- swapAux
- }
- if(integrate(Vectorize(f), lower = lowerBound, upper = upperBound)$value < 0)
- return(FALSE)
- }
- }
- return(TRUE)
- }
- else
- return(FALSE)
- },
- error = function(err){
- message(paste("eroare a functiei isPdf: ", err))
- return(FALSE)
- })
- }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement