[ADMB Users] error when imposing penalty on objective function
Jim Ianelli
Jim.Ianelli at noaa.gov
Thu Oct 7 08:31:29 PDT 2010
Hi,
In ADMB the indices of vector operations must be consistent. I.e., if you
subtract a (sub)vector of X which has dimension 1,10 e.g., :
X(2,10) - X(1,9);
will return an error. A shift operator is required (page 15-10 in my
manual). A shorthand way to do this is:
X(2,10) - ++X(1,9);
A long hand way is:
X(2,10) - X(1,9).shift(1);
Cheers,
Jim
James Ianelli
<http://www.afsc.noaa.gov/REFM/stocks/assessments.htm> REFM Division,
Alaska Fisheries Science Center
National Marine Fisheries Service, National Oceanic and Atmospheric
Administration
7600 Sand Point Way NE
Seattle WA 98115
Ph. 206 526 6510
From: users-bounces at admb-project.org [mailto:users-bounces at admb-project.org]
On Behalf Of Luis Ridao
Sent: Thursday, October 07, 2010 7:45 AM
To: users at admb-project.org
Subject: [ADMB Users] error when imposing penalty on objective function
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/5a95bb02/attachment.html>
More information about the Users
mailing list