[ADMB Users] Random effects model with integration in the likelihood

[Jeff Laake]

I'm attempting to incorporate a random effect into the scale of the
detection function for distance sampling which could be very useful for
fitting this type of data.  With Hans help I was successful when the
integral for the detection function could be solved analytically. However,
I'm having a problem incorporating adromb for numerical integration.  The
issue is that the function being integrated depends on the parameter sigma
which includes a random effect. adromb doesn't have an argument for
parameters of the function being integrated (fct below). Without random
effects, sigma was in the parameter section and was available globally but
that doesn't seem to work with random effects.  Below is my attempt at a
solution but it fails to compile when it hits the assignment for sigma. Any
suggestions/help would be much appreciated. You can find a pdf describing
what I'm trying to do at:
https://github.com/downloads/jlaake/ADMButils/distance_random_effect.pdf

DATA_SECTION
   init_int n;                        // number of distances
   init_number width;
   init_vector xs(1,n);               // distances
PARAMETER_SECTION
   init_number beta;                  // beta parameter for log-sigma;
   init_bounded_number sigeps(0.000001,5);
                                      // sigma for random effect;
   random_effects_vector u(1,n);      // random effect for scale
   objective_function_value f;        // negative log-likelihood
GLOBALS_SECTION
  #include <admodel.h>
  dvariable sigma;

PROCEDURE_SECTION
   int j;
// loop over each observation computing sum of log-likelihood values
// integral is sqrt(pi/2)*sigma; I dropped constant
   f=0;
   for (j=1;j<=n;j++)
   {
      ll_j(j,beta,sigeps,u(j));
   }

SEPARABLE_FUNCTION void ll_j(const int j, const dvariable& beta,const
dvariable& sigeps,const dvariable& u)
   dvariable eps=u*sigeps;
   int k;
   dvariable mu;
   sigma=exp(beta+eps);
   mu=adromb(&model_parameters::fct,0,width,8);
   f -= -0.5*square(u);
   f -= -log(mu) - 0.5*square(xs(j)/sigma);

FUNCTION dvariable fct(const dvariable& x)
// x is integration variable
// ifct is index for function read from data
   dvariable tmp;
   tmp=mfexp(-.5*x*x/square(sigma));
   return tmp;
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.admb-project.org/pipermail/users/attachments/20121109/f450f54a/attachment.html>