Сокращение вложенных циклов внутри функций с использованием lapply / sapply с несколькими параметрами и индексами - PullRequest
Новая
       27

Сокращение вложенных циклов внутри функций с использованием lapply / sapply с несколькими параметрами и индексами

0 голосов
/ 25 февраля 2019

У меня есть два вектора наблюдений, из которых я хочу вычислить, какие параметры лучше всего соответствуют этим наблюдениям.Поэтому я определил функцию правдоподобия.Я изо всех сил пытаюсь перейти от вложенных циклов к использованию Саппли.Моя цель - вернуть (50x1) вектор вероятностей наблюдения "бутона". 

> 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

Текущий код (который работает, хотя и медленно): 

# 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

Этот код дает(50x1) вектор с вероятностями наблюдения значений «бутонов», чего я и хочу.

Вот моя попытка перейти от цикла for к подходу sapply 

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)

Ошибка в pWEI2 (bud [i], theta1, theta2): отсутствует аргумент "theta1" без значения по умолчанию

Кажется, что аргументы не проходят весь путь через мой первыйфункция в моем коде.Будем благодарны за любые советы, как сделать так, чтобы мои аргументы передавались для всех функций!

1 Ответ

0 голосов
/ 28 февраля 2019

Мне удалось решить мою проблему. 

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))
...