[ADMB Users] error when imposing penalty on objective function

[Luis Ridao]

ADMB-help,

This is my first e-mail to this mailing-list. I'm new to ADMB and
trying to learn little by little.

I'm working on a catch-at-age model provided in the admb manual.
The first step was to import my personal stock data set. As a second step
I wish to impose a penalty on the objective function based on the
differences in the fishing mortality
coefficient ("log_fy_coff")

ADMB fails to run the model with the following message:

"Incompatible bounds in prevariable operator * (_CONST dvar_vector&
v1,_CONST dvar_vector& v2)"
The line causing the error is the  FUNCTION
"evaluate_the_objective_function" section at the of the file

# originally (from the manual admb.pdf)
FUNCTION evaluate_the_objective_function
  // penalty functions to ``regularize '' the solution
  f+=.01*norm2(log_relpop);
  avg_F=sum(F)/double(size_count(F));
  if (last_phase())
  {
    // a very small penalty on the average fishing mortality
    f+= .001*square(log(avg_F/.2));
  }
  else
  {
    f+= 1000.*square(log(avg_F/.2));
  }
  f+=0.5*double(size_count(C)+size_count(log_fy_coff))
    * log( sum(elem_div(square(C-obs_catch_at_age),.01+C))
    + 100*norm2(log_fy_coff)));     // ORIGINAL LINE

# the modification of the original code (last line)
FUNCTION evaluate_the_objective_function
  // penalty functions to ``regularize '' the solution
  f+=.01*norm2(log_relpop);
  avg_F=sum(F)/double(size_count(F));
  if (last_phase())
  {
    // a very small penalty on the average fishing mortality
    f+= .001*square(log(avg_F/.2));
  }
  else
  {
    f+= 1000.*square(log(avg_F/.2));
  }
  f+=0.5*double(size_count(C)+size_count(log_fy_coff))
    * log( sum(elem_div(square(C-obs_catch_at_age),.01+C))
    + 100*norm2(log_fy_coff(2,nyrs)-log_fy_coff(1,nyrs-1))); // MODIFIED
LINE

The .tpl file is provided below.
I'm running on a Windows machine.

########################################################################################
DATA_SECTION
  init_int nyrs                                 // the number of years odf
data
  init_int nages                                // the number of age
classess in the population
  init_matrix obs_catch_at_age(1,nyrs,1,nages)  // observed catch-at-age
data
  init_number M                                 // estimate of natural
mortality rate
  init_vector relwt(2,nages);                   // need to have relative
weight-at-age to calculate B2+
  vector ages(1,nages);                         // ages of data
  vector ages4plus(1,nages-1);                  // non-recruiting ages in
the population
  vector years(1,nyrs);                         // years of data
  int pred_year;                                // prediction year
INITIALIZATION_SECTION
  //log_q -1                                    // original -1
  log_popscale 5                                // original 5
PARAMETER_SECTION
  //init_number log_q(1)                        // log-catchability
  init_number log_popscale(1)                   // overall population
scaling parameter
  init_bounded_dev_vector log_sel_coff(1,nages-1,-15.,15.,2)        //
original log_sel_coff(1,nages-1,-15.,15.,2)
  init_bounded_dev_vector log_relpop(1,nyrs+nages-1,-15.,15.,2)     //
original log_relpop(1,nyrs+nages-1,-15.,15.,2)
  init_bounded_dev_vector log_fy_coff(1,nyrs,-.3,.3,3)              //
original log_fy_coff(1,nyrs,-2.,2.,3)
  vector log_sel(1,nages)
  vector log_fy(1,nyrs)
  vector log_initpop(1,nyrs+nages-1);
  matrix F(1,nyrs,1,nages)                     // instantaneous fishing
mortality
  matrix Z(1,nyrs,1,nages)                     // instantaneous total
mortality
  matrix S(1,nyrs,1,nages)                     // survival rate
  matrix N(1,nyrs,1,nages)                     // predicted numbers-at-age
  matrix C(1,nyrs,1,nages)                     // predicted catch-at-age
  objective_function_value f
  number recsum
  number initsum
  sdreport_number avg_F
  sdreport_vector predicted_N(2,nages)
  sdreport_vector ratio_N(2,nages)
  // changed from the manual because adjusted likelihood routine doesn't
  // work
  likeprof_number pred_B

PRELIMINARY_CALCS_SECTION
  ages.fill_seqadd(3,1);           // vector of ages
  ages4plus.fill_seqadd(4,1);      // vector of non recruiting ages
  years.fill_seqadd(1975,1);       // fill vector of years with years
  pred_year=years[nyrs]+1;         // year of prediction = last_year +1
PROCEDURE_SECTION
  // example of using FUNCTION to structure the procedure section
  get_mortality_and_survivial_rates();

  get_numbers_at_age();

  get_catch_at_age();

  evaluate_the_objective_function();

FUNCTION get_mortality_and_survivial_rates
  int i, j;
   //calculate the selectivity from the sel_coffs ---------------------
  for (j=1;j<nages;j++)
  {
  log_sel(j)=log_sel_coff(j);
  }
  //the selectivity is the same for the last two age classes
  log_sel(nages)=log_sel_coff(nages-1);

  for (i=1;i<=nyrs;i++)
  {
    log_fy(i)=log_fy_coff(i);
  }
  F=outer_prod(mfexp(log_fy),mfexp(log_sel));  //
F=outer_prod(mfexp(log_fy),mfexp(log_sel)) ;
F=outer_prod(mfexp(log_q)*effort,mfexp(log_sel));
  Z=F+M;
  // get the survival rate
  S=mfexp(-1.0*Z);

FUNCTION get_numbers_at_age
  int i, j;

  log_initpop=log_relpop+log_popscale;

  for (i=1;i<=nyrs;i++)
  {
    N(i,1)=mfexp(log_initpop(i));
  }
  for (j=2;j<=nages;j++)
  {
    N(1,j)=mfexp(log_initpop(nyrs+j-1));
  }
  for (i=1;i<nyrs;i++)
  {
    for (j=1;j<nages;j++)
    {
      N(i+1,j+1)=N(i,j)*S(i,j);
    }
  }
  // calculated predicted numbers at age for next year
  for (j=1;j<nages;j++)
  {
    predicted_N(j+1)=N(nyrs,j)*S(nyrs,j);
    ratio_N(j+1)=predicted_N(j+1)/N(1,j+1);
  }
  // calculated predicted Biomass for next year for
  // adjusted profile likelihood
  pred_B=(predicted_N * relwt);

FUNCTION get_catch_at_age
  C=elem_prod(elem_div(F,Z),elem_prod(1.-S,N));

FUNCTION evaluate_the_objective_function
  // penalty functions to ``regularize '' the solution
  f+=.01*norm2(log_relpop);
  avg_F=sum(F)/double(size_count(F));
  if (last_phase())
  {
    // a very small penalty on the average fishing mortality
    f+= .001*square(log(avg_F/.2));
  }
  else
  {
    f+= 1000.*square(log(avg_F/.2));
  }
  f+=0.5*double(size_count(C)+size_count(log_fy_coff))
    * log( sum(elem_div(square(C-obs_catch_at_age),.01+C))
    + 100*norm2(log_fy_coff(2,nyrs)-log_fy_coff(1,nyrs-1)));
########################################################################################
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.admb-project.org/pipermail/users/attachments/20101007/e1d0d854/attachment.html>