2016-09-28 2 views
0

J'essaie d'estimer les trois paramètres a, b0 et b1 avec la fonction optim(). Mais j'obtiens toujours l'erreur: Erreur dans optim (par = c (1, 1, 1), fn = logweibull, méthode = "L-BFGS-B",: L-BFGS-B a besoin de valeurs finies de 'fn 'R optim() L-BFGS-B a besoin de valeurs finies de 'fn' - Weibull

t<-c(6,6,6,6,7,9,10,10,11,13,16,17,19,20,22,23,25,32,32,34,35,1,1,2,2,3,4,4,5,5,8,8,8,8,11,11,12,12,15,17,22,23) 
d<-c(0,1,1,1,1,0,0,1,0,1,1,0,0,0,1,1,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1) 
X<-c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0) 

logweibull <- function (a,b0,b1) {a <- v[1];b0 <- v[2]; b1 <- v[3]; 
sum (d*log(t^a*exp(b0+X*b1)-t^a*exp(b0+X*b1))) + sum (d + log((a*t^(a-1))/t^a)) } 

v<-c(1,1,1) 

optim(par=c(1,1,1) ,fn = logweibull, method = "L-BFGS-B",lower = c(0.1, 0.1,0.1), upper = c(100, 100,100),control = list(fnscale = -1)) 

pouvez-vous me aider? savez-vous ce que je fait de mal?

+0

'logweibull <- function (a, b0, b1) {a <- v [1]; b0 <- v [2]; b1 <- v [3]; sum (d * log (t^a * exp (b0 + X * b1)) - t^a * exp (b0 + X * b1)) + somme (d * log ((a * t^(a-1))/t^a))} 'Merci, il y a eu une erreur et il manque des parenthèses. Mais maintenant, la fonction ne fonctionne toujours pas correctement. Je vais récupérer seulement les valeurs de départ en tant que paramètres estimés. Avez-vous une idée de ce qui ne va pas? – Hans

Répondre

2

vous pouvez également envisager

(1) passer les variables de données supplémentaires à la fonction objective ainsi que les paramètres

(2) passer la fonction de gradient (ajouté t il gradient fonction)

(3) la fonction objectif original peut être encore simplifiée (comme ci-dessous)

logweibull <- function (v,t,d,X) { 
    a <- v[1] 
    b0 <- v[2] 
    b1 <- v[3] 
    sum(d*(1+a*log(t)+b0+X*b1) - t^a*exp(b0+X*b1) + log(a/t)) # simplified function 
} 

grad.logweibull <- function (v,t,d,X) { 
    a <- v[1] 
    b0 <- v[2] 
    b1 <- v[3] 
    c(sum(d*log(t) - t^a*log(t)*exp(b0+X*b1) + 1/a), 
    sum(d-t^a*exp(b0+X*b1)), 
    sum(d*X - t^a*X*exp(b0+X*b1))) 
} 

optim(par=c(1,1,1), fn = logweibull, gr = grad.logweibull, 
     method = "L-BFGS-B", 
     lower = c(0.1, 0.1,0.1), 
     upper = c(100, 100,100), 
     control = list(fnscale = -1), 
     t=t, d=d, X=X) 

avec sortie

$par 
[1] 0.2604334 0.1000000 0.1000000 

$value 
[1] -191.5938 

$counts 
function gradient 
     10  10 

$convergence 
[1] 0 

$message 
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH" 

En outre, ci-dessous est une comparaison entre la convergence avec et sans fonction de gradient (avec différence finie). Avec une fonction de gradient explicite, il faut 9 itérations pour converger vers la solution, alors que sans elle (avec différence finie), il faut 126 itérations pour converger.

enter image description here

+0

Mais idéalement, l'espérance conditionnelle ne devrait-elle pas être calculée en utilisant le pdf conditionnel? jetez un oeil à ceci: http://stats.stackexchange.com/questions/12843/generating-random-samples-from-a-custom-distribution, comment ils échantillonnent à partir d'un pdf arbitraire (en calculant la cdf puis en trouvant racines), l'échantillonnage de Monte Carlo ne devrait-il pas être fait à partir de ce pdf tronqué? –