[R] maxLik pakage

[spencerg]

  1. Have you worked through the examples in the maxLik help page? Your 
example is sufficiently complicated that I hesitate to try it myself, 
especially since I see characters in your email that are not simple 
ASCII. If you can get the examples in the maxLik help page to work, 
identify the differences between your example and those on the help 
page, then generate intermediate examples to test your understanding 
until you either solve your problem or you find a simpler example that 
you don't understand and then post that to this list.


2. Have you checked your gradient with the "compareDerivatives" function 
in the maxLik package? If your gradient function is not correct, this 
should help you understand and hopefully fix the problem.


Hope this helps.
Spencer

amene kheradmandi wrote:
> Hi all;
> I recently have been used 'maxLik' function for maximizing G2StNV178 function with gradient function gradlik; for receiving this goal, I write the following program; but I have been seen an error  in calling gradient  function;
> The maxLik function can't enter gradlik function (definition of gradient function); I guess my mistake is in line ******** ,that the vector  ‘h’ is defined.
> Â 
> I need your advices for resolve the error;
> Sincerely yours
> A. Kheradmandi
> ##########data vector
> x=c(0,0,14,7,0,22,2,18,10,29,
> 16,11,59,6,33,5,1,0,28,17,
> 28,19,10,43,9,17,42,2,22,6,
> 14,0,59,9,11,5,11,37,15,26,
> 57,8,16,22,5,13,19,12,34,8,28,35,13,11,59,
> 8,7,20,36,39,14,31,18,32,30,
> 29,9,2,5,3,19,27,76,4,32,
> 24,13,3,13,24,41,20,0,13,26,
> 0,23,19,10,14,14,13,15,18,29,
> 29,0,20,0,22,13,17,12,29,0,
> 0,16,55,11,19,9,2,9,32,16,
> 55,10,59,8,65,39,1,8,2,7,
> 37,8,0,5,8,0,7,0,20,10,
> 11,3,0,7,31,14,18,3,4,34,
> 32,4,9,9,8,36,44,0,9,27,28,55,72,12,1,
> 9,0,32,0,0,2,15,5,6,17,
> 63,61,9,15,15,0,2)
> #########goal is  found unique maximum for  5-parameter function with using  "maxLik" function
> #####the function in latex commands is as following:
> Â \begin{eqnarray*}
> &&\ell(\xi,\omega,\nu,\lambda_1,\lambda_2)=n\log2-n\log\omega+n\log\Gamma(\frac{\nu+1}{2})-\frac{n}{2}\log(\nu\pi)\\
> &-&n\log\Gamma(\frac{\nu}{2})-\frac{\nu+1}{2}\sum_{i=1}^{n}\log(1+\frac{(x_i-\xi)^2}{\omega^2\nu})+\sum_{i=1}^{n}\log\Phi(\lambda_1\frac{(x_i-\xi)}{\sqrt{\omega^2+\lambda_2(x_i-\xi)^2}})
> \end{eqnarray*}
> ##############
> G2StNV178<-function(a){
> require(maxLik)
> II=0
> nu<-(a[1])
> lambda1<-a[2]
> lambda2<-a[3]
> ksi<-a[4]
> omega<-a[5]
> II<-log(prod((2*dt((x-ksi)/omega,nu)*pnorm((lambda1*(x-ksi)/omega)/sqrt(1+lambda2*((x-ksi)/omega)^2)))/omega))
> II
> }
> ########definition of gradient  function
> gradlik<- function(a){
> nu<-a[1]
> lambda1<-a[2]
> lambda2<-a[3]
> ksi<-a[4]
> omega<-a[5]Â  
> #*******#
> h<-c(f1=n*digamma((nu+1)/2)-n/(2*nu)-n*digamma(nu/2)-1/2*sum(log(1+(x-ksi)^2/(omega^2nu)))+(nu+1)/2*sum((x-ksi)^2/(omege^2*nu^2+nu*(x-ksi)^2)),
> f2=sum((x-ksi)*dnorm(lambda1*(x-ksi)/sqrt(omega^2+lambda2*(x-ksi)^2))/(sqrt(omega^2+lambda2*(x-ksi)^2)*pnorm(lambda1*(x-ksi)/sqrt(omega^2+lambda2*(x-ksi))))), 
> f3=sum((lambda1*(x-ksi)^3*dnorm(lambda1*(x-ksi)/sqrt(omega^2+lambda2*(x-ksi)^2)))/( 2*omega^3*pnorm(lambda1*(x-ksi)/sqrt(omega^2+lambda2*(x-ksi)^2))
> *sqrt((1+lambda2*(x-ksi)^2/omega^2 )^3))),
> f4=(nu+1)*sum((x-ksi)/(omega^2*nu+(x-ksi)^2))-sum((lambda1*dnorm((lambda1*(x-ksi))/(omega^2+lambda2*(x-ksi)^2)))/(omega*sqrt(1+lambda2*(x-ksi)^2/(omega^2))
> *pnorm(lambda1*(x-ksi)/(sqrt(omega^2+lambda2*(x-ksi)^2))))), 
> f5=-n/omega+(nu+1)*n*sum((x-ksi)^2/(nu*omega^3+omega*(x-ksi)^2))-sum((lambda1*(x-ksi)*omega*dnorm(lambda1*(x-ksi)/(sqrt(omega^2+lambda2*(x-ksi)^2))
> /(pnorm(lambda1*(x-ksi)/(sqrt(omega^2+lambda2*(x-ksi)^2))*sqrt((omega^2+lambda2*(x-ksi)^2)^3)))))))
> }
> n<- length(x) # the final command is :
> res<- maxLik(G2StNV178,grad=gradlik,start=c(4.666484 ,4603.399055 , 4603.399055Â Â  ,-0.016674 ,19.055865),tol=1e-16,iterlim =3000)
> summary(res)
>
>
>       
> 	[[alternative HTML version deleted]]
>
>   
> ------------------------------------------------------------------------
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>