[ADMB Users] lgam in 10.0
dave fournier
davef at otter-rsch.com
Sun May 15 11:42:40 PDT 2011
It raining again (still).
Thought I would look at this.
Below is the whole TPL file as I have modified it to get it to run.
first problme is that posfun must have a non negative second argument.
This is where the "positiveness" kicks in. If you make it 0.0 the
curvature is infinite. Best is to tune it as the model gets fitted.
Call it eps2 and begin with a bigger number in the first phase.
This is a typical technique to make the objective function well
behaved until you get near the solution.
x1=posfun(prevhr-starthr,eps2,fpen);
x2=posfun(hr-starthr,eps2,fpen);
Then
f=sum(elem_prod(eclosed,log(eps1+p))+elem_prod((avail-eclosed),
log(eps1+(1-p))));
will cause you trouble if p is too close to 0 or 1. so add a tiny number
and tune it as well.
Note the order of addition (eps1+(1-p)) ans figure out why (eps1+1-p)
may not work.
Then add a kludge variable with the number of phases that you want. I
always name if
kkludge so that I l know what its purpose is.
> f+=square(kkludge);
Now run the model and tune the values of eps1 eps2 to be small but the
model still converges.
I got into this and ended up with 5 phases. the value of eps1 has a
large effect.
I ran with the CLO of -crit 1.e-8.
DATA_SECTION
init_int nobs
init_int nttt
init_vector eclosed(1,nobs)
init_vector avail(1,nobs)
init_vector hr(1,nobs)
init_vector prevhr(1,nobs)
init_matrix ttt(1,nobs,1,nttt)
PARAMETER_SECTION
init_vector tttcoefs(1,nttt)
init_number k
init_number starthr
init_number kkludge(5);
vector loglam(1,nobs)
vector p(1,nobs)
vector lambda(1,nobs)
vector x1(1,nobs)
vector x2(1,nobs)
objective_function_value f
PROCEDURE_SECTION
double eps1=1.e-13;
double eps2=0.001;
if (current_phase()<2)
{
eps1=1.e-2;
eps2=.02;;
}
else if (current_phase()<3)
{
eps1=1.e-8;
}
else if (current_phase()<4)
{
eps1=1.e-10;
}
else if (current_phase()<5)
{
eps1=1.e-12;
}
loglam=ttt*tttcoefs;
lambda=exp(loglam);
dvariable fpen=0.0;
x1=posfun(prevhr-starthr,eps2,fpen);
x2=posfun(hr-starthr,eps2,fpen);
p=exp(-pow(elem_div(x1,lambda),k))-exp(-pow(elem_div(x2,lambda),k));
f=sum(elem_prod(eclosed,log(eps1+p))+elem_prod((avail-eclosed),
log(eps1+(1-p))));
f+=square(kkludge);
TOP_OF_MAIN_SECTION
arrmblsize = 2147483647;
More information about the Users
mailing list