[ADMB Users] error when imposing penalty on objective function
Steve Martell
s.martell at fisheries.ubc.ca
Thu Oct 7 08:47:39 PDT 2010
Luis,
The error is occurring because you are trying to subtract a vector from another vector where the bounds are different. (log_fy_coff(2,nyrs)-log_fy_coff(1,nyrs-1)) This is not allowed. You can use the .shift function to change the array bounds (see the ADMB manual at the bottom of page 168, or 15-10 for more information). Also it looks like your trying to minimize the first difference in the log_fy_coff, there is an admb function called "first_difference" that will accomplish what you are trying to do. e.g. norm2(first_difference(log_fy_coff)).
Steve
On 2010-10-07, at 7:44 AM, Luis Ridao wrote:
> 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)));
> ########################################################################################
> _______________________________________________
> Users mailing list
> Users at admb-project.org
> http://lists.admb-project.org/mailman/listinfo/users
Steve Martell
s.martell at fisheries.ubc.ca
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.admb-project.org/pipermail/users/attachments/20101007/55ff2560/attachment.html>
More information about the Users
mailing list