[R] declaring negative log likelihood of a distribution

chamilka cmanoj4 at gmail.com
Wed Jul 11 15:11:07 CEST 2012


Hi everyone!
I already posted 
http://r.789695.n4.nabble.com/Declaring-a-density-function-with-for-loop-td4635699.html
a question  on finding density values of a new Binomial like distribution
which has the following pmf:
http://r.789695.n4.nabble.com/file/n4636134/kb.png 
 Thank fully
http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=user_nodes&user=124474 
Berend Hasselman  and
http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=user_nodes&user=125777 
David Carlson  cleared out my problems in a timely manner.
 
 Now I have a question on estimating the parameters of this distribution
when a data set is provided. 
 For example, I have a data from a journal paper:

*> values<- c(0,1,2,3,4,5,6,7) #Binomial variable values
> frequency<-c(47,54,43,40,40,41,39,95) #frequency counts of the  above
> binomial values, totally I have 399 samples
 *

 Now to estimate the parameters, I first declared the negative log
likelihood function of this new distribution as below(I used two methods,
one with a loop and one without the loop)

#without any loops
*> Loglik.newdis1<- function(x,a,b,fre,n){
      KBLL1<-sum(fre*(log(a*b)+lchoose(n,x)+log(sum((-1)**(0:1000) *
choose(b-1, 0:1000) * beta(x+a+a*0:1000, n-x+1)))))
      return(-KBLL1)
      }*
     #here a and b  are the parameters to be estimated and x=binomial
values, fre=frequencies and n=binomial trials
     #and since the inner series is a convergent infinite series, a large
number like 1000 is defined as the maximum
     
 #with for loop   
*>Loglik.newdis2<-function(x,a,b,fre,n,imax) { 
      term<-0 
      for (i in 0:imax) { 
           term=term+(((-1)**i)*(choose(b-1,i))*(beta(x+a+a*i,n-x+1))) 
           }
      dens=a*b*choose(n,x)*term 
      KBLL2<-sum(fre*log(dens))
      return(-KBLL2) 
      }  *
      
 #Now to estimate the parameters, I used 
http://rgm2.lab.nig.ac.jp/RGM2/func.php?rd_id=bbmle:mle2 mle2(from the
package bbmle) , which is a wrapper around the optim function

*>estimates1=mle2(Loglik.newdis1, start=list(a=1,b=1),
data=list(x=values,fre=frequency, n=7, imax=2000)) *
  Error in optim(par = c(1, 1), fn = function (p)  : 
  non-finite finite-difference value [1]
In addition: There were 50 or more warnings (use warnings() to see the first
50)
/#this gives the above error message/

>estimates2=mle2(Loglik.newdis2, start=list(a=1,b=1),
data=list(x=values,fre=frequency, n=7, imax=2000))
 > estimates2

Call:
mle2(minuslogl = Loglik.newdis2, start = list(a = 1, b = 1), 
    data = list(x = values, fre = frequency, n = 7, imax = 2000))

Coefficients:
        a         b 
0.7574478 0.6399205 

Log-likelihood: -815.48 

     /#here the values are not equal to the vallues provided in the journal
paper
     #also thee inner series doesn't converge even for very large number
     #values given in the paper are 
     a= 0.70 and b= 0.59 and -loglikelihood = 813.6939/
*
Please point out my mistakes and help me on this r-code?
Thank you very much!!*

--
View this message in context: http://r.789695.n4.nabble.com/declaring-negative-log-likelihood-of-a-distribution-tp4636134.html
Sent from the R help mailing list archive at Nabble.com.



More information about the R-help mailing list