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

John Sibert sibert at hawaii.edu
Fri Nov 9 11:11:28 PST 2012


Try getting rid of the GLOBALS_SECTION and declaring sigma as a number 
in the PARAMETER_SECTION

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
    number 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 tm

John Sibert
Emeritus Researcher, SOEST
University of Hawaii at Manoa

Visit the ADMB project http://admb-project.org/

On 11/09/2012 07:00 AM, Jeff Laake wrote:
> 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;
>
>
>
>
> _______________________________________________
> Users mailing list
> Users at admb-project.org
> http://lists.admb-project.org/mailman/listinfo/users




More information about the Users mailing list