Skip to content

help with the maxBHHH routine

6 messages · Andrew Robinson, Ravi Varadhan, Arne Henningsen +1 more

Original
Andrew Robinson
# 
I suggest that you provide some commented, minimal, self-contained,
reproducible code.

Cheers

Andrew
On Wed, May 04, 2011 at 02:23:29AM +0530, Rohit Pandey wrote:

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Ravi Varadhan
    

    

      
      #
    

  


  

      
maxBHHH is *not* an in-built R function.  It is in a distributed package called "maxLik".  Always tell us which package is being used so that it is easier for us to help you.

The error message says that the gradient function is returning a 10 x 2 matrix, whereas you say that you have 1000's of observations and 821 parameters.  Show us a simplified version of your problem.  It is difficult to help you without seeing an example. 

Also, try running w/o the gradient and see how it works.  This is not the answer you wanted, but we cannot help you w/o seeing your example.

Ravi.

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Arne Henningsen
    

    

      
      #
    

  


  

      
Dear Rohit

        
On 3 May 2011 22:53, Rohit Pandey <rohitpandey576 at gmail.com> wrote:

            

      
      
      
    

        
Please make yourself familiar with the BHHH algorithm and read the
documentation of maxBHHH: it says about argument "grad":

"[...] If the BHHH method is used, ?grad? must return a matrix, where
rows correspond to the gradient vectors of individual observations and
the columns to the individual parameters.[...]"

More information of the maxLik package is available at:
http://dx.doi.org/10.1007/s00180-010-0217-1

Best regards,
Arne

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Rohit Pandey
    

    

      
      #
    

  


  

      
Hi Andrew, Ravi and Arne,

Thank you so much for your prompt replies. I see that all of you mention the
need for simple, reproducible code. I had thought of doing this, but the
functions I was using for the observation level gradient and likelihood
function were very long. I will paste them below here.

Also, sorry for the ambiguity with the "1000's of observations and 821
parameters" on the one hand and the 10 * 2 matrix on the other. The latter
is a toy data set and the former is the real data set I ultimately hope to
apply this routine to once it works. Also, sorry for not mentioning the fact
that the maxBHHH function I am using is from the maxLik package (thanks,
Ravi for pointing out).

So, the code that is giving me the errors is:

maxBHHH(logLikALS4,grad=nuGradientC4,finalHessian="BHHH",start=prm,iterlim=2)

and

maxBHHH(logLikALS4,grad=nuGradientC4,finalHessian="BHHH",start=prm,iterlim=2)
Where nuGradientC4 returns a 2*10 matrix and nuGradientC5 a 10*2 matrix
(there are 10 parameters and 2 observations).

I have attached the required functions in the .R file.

These make for some pretty long code, but all you have to do is either load
the file or paste the contents into your R console (and maybe see that
they're returning what they're supposed to). I'm sorry I couldn't think of a
way to come up with a shorter version of this code (I tried my best).

Once you load the file, you should see the following:


#The observation level likelihood function

            

      
      
      
    

        
1          2
-0.6931472 -0.6931472

#The observation level gradients

            

      
      
      
    

        
1          2         3 4          5          6
7          8         9        10
2 -0.3518519  0.3518519 0.0000000 0 -0.1481481 -0.1666667  0.1481481
0.1666667 0.0000000 0.0000000
4  0.0000000 -0.3518519 0.3518519 0  0.0000000  0.0000000 -0.1666667
-0.1481481 0.1666667 0.1481481
Warning messages:
1: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'
2: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'

            

      
      
      
    

        
2          4
1  -0.3518519  0.0000000
2   0.3518519 -0.3518519
3   0.0000000  0.3518519
4   0.0000000  0.0000000
5  -0.1481481  0.0000000
6  -0.1666667  0.0000000
7   0.1481481 -0.1666667
8   0.1666667 -0.1481481
9   0.0000000  0.1666667
10  0.0000000  0.1481481
Warning messages:
1: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'
2: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'

Ignore the warning messages.

The errors are:

            

      
      
      
    

        
maxBHHH(logLikALS4,grad=nuGradientC4,finalHessian="BHHH",start=prm,iterlim=2)
Error in checkBhhhGrad(g = gr, theta = theta, analytic = (!is.null(attr(f,
:
  the matrix returned by the gradient function (argument 'grad') must have
at least as many rows as the number of parameters (10), where each row must
correspond to the gradients of the log-likelihood function of an individual
(independent) observation:
 currently, there are (is) 10 parameter(s) but the gradient matrix has only
2 row(s)
In addition: Warning messages:
1: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'
2: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'

 and:

            

      
      
      
    

        
maxBHHH(logLikALS4,grad=nuGradientC5,finalHessian="BHHH",start=prm,iterlim=2)
Error in gr[, fixed] <- NA : (subscript) logical subscript too long
In addition: Warning messages:
1: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'
2: In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'

Again, thanks for your patience and help.

Rohit


On Wed, May 4, 2011 at 4:44 AM, Andrew Robinson <

        
A.Robinson at ms.unimelb.edu.au> wrote:

        

            

      
      
      
    

        

      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Andrew Robinson
    

    

      
      #
    

  


  

      
Hi Rohit,

actually, the request for simple reproducible code means that you have
to find the simplest possible representation of the problem.  

What happens if you simplify the observation level gradient and the
likelihood function?  Eg to trivial examples?  If you still get the
error, then simplify it futher.  If you get the error with the
simplest possible problem, then share it.  If you don't , then try to
figure out what the changes were that resolved the problem, and scale
those back up to your original problem.

Does that make sense?

Cheers

Andrew

        
On Thu, May 05, 2011 at 03:22:55AM +0530, Rohit Pandey wrote: