Reducing nested loops within functions using lapply/sapply with several parameters and indexes

Question

I have two vectors of observations, of which I want to calculate what parameters which best fit these observations. Therefore, I have defined a likelihood function. I am struggling to go from nested for loops to using sapply. My goal is to return a (50x1) vector of probabilities of observing "bud".

> bud
 [1] 1.72 1.12 1.01 1.14 1.56 1.04 1.00 1.39 1.38
[10] 0.99 1.66 1.55 2.26 0.86 1.45 0.90 1.58 1.09
[19] 0.97 0.98 1.16 1.20 1.17 1.59 1.00 0.92 1.24
[28] 0.98 1.29 1.36 1.80 1.42 1.67 1.13 0.76 0.93
[37] 1.03 1.29 1.09 1.64 1.27 1.14 1.60 1.00 0.91
[46] 1.40 1.15 2.03 1.34 1.37
> deltakere
 [1] 8 5 7 6 9 5 5 9 9 6 5 8 8 6 6 5 7 6 7 7 9 6 6 8
[25] 8 6 6 5 9 9 6 7 9 5 5 6 6 8 7 9 8 6 9 8 5 8 8 6
[49] 9 7

The current code (which works, albeit slowly):

# first equation in my likelihood function
eqOne <- function(i,theta1, theta2, lambda){
        vektorN = (1:100)
        loopVerdi = 0
            if(deltakere[i]>=2){        
                for (N in deltakere[i]:20000) {
                    density = ((N)*(N-1)*(pWEI2(bud[i], theta1, theta2))^(N-2) *
                        pWEI2(bud[i], theta1, theta2, lower.tail = FALSE) *
                        dWEI2(bud[i], theta1, theta2)) /
                        (1-(pWEI2(r, theta1, theta2))^N)
                    ngittN =
                        dbinom(deltakere[i], size = N, 
                               prob = pWEI2(r, theta1, theta2, lower.tail = FALSE))
                    sN = 
                        dpois(N, lambda)
                    loopVerdi = loopVerdi + (density * ngittN * sN)}
            return(loopVerdi)}}

# here is the second equation in the likelihood estimate.
eqTwo <- function(theta1, theta2, lambda) {
        firstPart <- 1:length(bud)
            for (i in 1:length(bud)){
                firstPart[i] <- eqOne(i,theta1, theta2, lambda)}
        return(firstPart)}

Vector <- eqTwo(1,2,7)
> Vector
 [1] 0.071126958 0.174198273 0.125381705 0.202950594
 [5] 0.116409283 0.156245797 0.143071864 0.156875579
 [9] 0.157670857 0.139149430 0.052085368 0.138248017
[13] 0.001884680 0.073326986 0.161989881 0.103352207
[17] 0.125124465 0.186201434 0.103156346 0.108636007
[21] 0.116075468 0.214339615 0.209921459 0.121433257
[25] 0.084230719 0.102485943 0.216127393 0.135710730
[29] 0.153834163 0.158552633 0.036673600 0.194744214
[33] 0.079498794 0.175511921 0.049722001 0.107646750
[37] 0.159465393 0.198467993 0.169062724 0.089375280
[41] 0.196173239 0.202950594 0.102923165 0.084230719
[45] 0.107521392 0.190136547 0.159243477 0.008274964
[49] 0.158457042 0.210361176

This code yields a (50x1) vector with probabilities of observing the "bud" values, which is what I want.

Here is my attempt on going from the for loop to a sapply approach

funDeltakere <- function(i, theta1, theta2, lambda){ 
    function(N, theta1, theta2, lambda){
                    density = ((N) *(N-1)*(pWEI2(bud[i], theta1, theta2))^(N-2) *
                        pWEI2(bud[i], theta1, theta2, lower.tail = FALSE) *
                        dWEI2(bud[i], theta1, theta2)) /
                        (1-(pWEI2(r, theta1, theta2))^N)
                    ngittN =
                        dbinom(deltakere[i], size = N, 
                               prob = pWEI2(r, theta1, theta2, lower.tail = FALSE))
                    sN = dpois(N, lambda)
                    return(density * ngittN * sN)}}

eqOne <- function(i,theta1, theta2, lambda){
            if(deltakere[i]>=2){        
                loopVerdi <- sapply(deltakere[i]:20000, funDeltakere(i, theta1, theta2, lambda))
            return(loopVerdi)}} 

eqTwo <- function(theta1, theta2, lambda) {
            firstPart <- 1:length(bud)
                for (i in 1:length(bud)){
                    firstPart[i] <- eqOne(i,theta1, theta2, lambda)
                }
            return(firstPart)
            }
Vector <- eqTwo(1,2,7)

Error in pWEI2(bud[i], theta1, theta2) : argument "theta1" is missing, with no default

It seems like the arguments are not passing all the way through my first function in my code. Any tips on how to make sure my arguments passes for all the functions would be appreciated!

Problems presented only with broken code are not particularly useful. You need to specify the problem in natural language, — IRTFM, Feb 25 '19 at 21:35
If you make this issue fully reproducible it will be much easier for us to help you. Check out this post for guidlines: https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example Also you can use the reprex-package to make sure your answer is reproducible: https://reprex.tidyverse.org/ — FilipW, Feb 26 '19 at 10:00
I have tried to make it reproducible and given explanation for what I want with the code. — Ola, Feb 26 '19 at 10:34

score 0 · Answer 1 · answered Feb 28 '19 at 09:30

I managed to solve my problem.

bud <- runif(10, 1.1,1.55)
deltakere <- round(runif(10,5,9))
r=0.5

foo_inner <- function(i, theta1, theta2, lambda){
    function(N){
    density = (
        (N) *
        (N-1) * 
        (pWEI2(bud[i], theta1, theta2))^(N-2) *
        pWEI2(bud[i], theta1, theta2, lower.tail = FALSE) *
        dWEI2(bud[i], theta1, theta2)) /
        (1-(pWEI2(r, theta1, theta2))^N)

    ngittN =
        dbinom(deltakere[i], size = N, 
               prob = pWEI2(r, theta1, theta2, lower.tail = FALSE))
    sN = 
        dpois(N, lambda)

    return(density * ngittN * sN)
    }}

foo_outer <- function(theta1, theta2, lambda){
            function(i){
            listeObs <- sapply(deltakere[i]:20000, foo_inner(i, theta1, theta2, lambda))
            return(sum(listeObs))
            }}

eqThree <- function(theta1, theta2, lambda){
    secondPart <- pbsapply(1:length(bud), foo_outer(theta1, theta2, lambda))

    LL <- -sum(log(secondPart))
    return(LL)
}

result_mle <- mle(minuslogl = eqThree, start=list(theta1 = 1,
                                                  theta2 = 2,
                                                  lambda = 7),
                  method="L-BFGS-B", lower=c(0.1,0.1,5),
                  nobs = length(bud))

Reducing nested loops within functions using lapply/sapply with several parameters and indexes

1 Answers1