[ADMB Users] next try: "Incompatible shapes in df1b2vector functionoperator +"
Bolker,Benjamin Michael
[Continuing with my saga of trying to get a
nonlinear size-density functional response model with
random effects running ...]
It does look like initial values for the random effects
vector are indeed required. Is this actually documented somewhere
that I missed ... ??
I got simple2.tpl to run OK, so I am not completely hopeless (maybe).
Re Dave Fournier's previous message: I agree that there are potential issues with the scale of
the random effects/ bounds on the probabilities, but I think that in reality for this problem
the random-effects variance is quite small, so I'm hoping to get away
(using good starting values and bounding the random effects variance)
with using this form for now. Or I will use posfun() and penalties to make things
behave. If it seems likely that this issue is causing my current problems I will
go ahead and try to address it now -- otherwise I would like to try to get
the model running with a highly constrained RE variance, then relax
the constraints a bit and see if I run into trouble.
I have instrumented my TPL file a little bit so that it tells
me what calculations it is doing. It makes a little bit of headway
(prints out finite values for the function value, claims to have
gotten to at least iteration 50 -- although if I'm reading the output
correctly some of the parameter values are negative despite being
bounded below at zero??? never mind ...)
It merrily prints out my debug statements for the calculations
a few hundred times, then gets to
inner maxg = 0.103391
Inner second time = 0 Inner f = 11.6341
Newton raphson 1
and then stops with
Incompatible shapes in df1b2vector functionoperator +
??? hard to tell where this is coming from ... ???
The attached data file is very small (only the first 10 records),
but I get similar results with the full data file.
I could try simulating well-behaved data, but (a) these data are
pretty well-behaved anyway and (b) I have a feeling that it won't
make a difference (but will give it a shot if someone thinks it will).
Any ideas?
sincerely
Ben Bolker
PIN file:
------------------
# "mccoypred5.pin" produced by pin_write() from ADMButils; Thu Aug 6 01:39:03 2009
# c
0.506685
# d
4.08488
# h
0.008644705
# g
3.019749
# sigma_c
0.1
# u
0 0 0 0 0 0
DAT file:
----------------
# "mccoypred5.dat" produced by dat_write() from ADMButils; Thu Aug 6 01:39:03 2009
# nobs
10
# nblock
6
# killed
3 0 1 3 1 0 1 2 2 2
# size
8.830333 9.155 7.991667 8.000833 8.140333 8.711167 16.5265 14.1508 18.8445 14.87525
# initial
6 6 6 6 6 6 6 6 6 6
# Z
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0 0
0 0 0 1 0 0
0 0 0 0 1 0
0 0 0 0 0 1
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0 0
0 0 0 1 0 0
TPL file:
DATA_SECTION
init_int nobs // # of observations
init_int nblock // # of blocks
init_vector killed(1,nobs) // # killed per trial
init_vector size(1,nobs) // size of individuals
init_vector initial(1,nobs) // starting density (# individuals)
init_matrix Z(1,nobs,1,nblock) // random-effects model matrix
PARAMETER_SECTION
init_bounded_number c(0.0,1.0,1)
init_bounded_number d(0.0,10,1)
init_bounded_number h(0.0,1.0,1)
init_bounded_number g(0.0,10,1)
init_bounded_number sigma_c(0.00001,0.1,1)
random_effects_vector u(1,nblock,2)
vector prob(1,nobs)
vector cvec(1,nobs)
vector a1vec(1,nobs)
objective_function_value f
PROCEDURE_SECTION
cout << "cvec calc\n";
cvec = sigma_c*(Z*u)+c;
a1vec = elem_prod(pow(size/d,g),exp(-size/d));
cout << "prob calc\n";
prob = 1/(1/(elem_prod(cvec,a1vec))+h*initial);
cout << "loglik calc\n";
f -= sum( log_comb(initial,killed)+
elem_prod(killed,log(prob))+
elem_prod(initial-killed,log(1-prob)));
GLOBALS_SECTION
#include <admodel.h>
#include <df1b2fun.h>
#include <adrndeff.h>
df1b2vector pow(const df1b2vector& v,const df1b2variable & x)
{
int mmin=v.indexmin();
int mmax=v.indexmax();
df1b2vector tmp(mmin,mmax);
for (int i=mmin;i<=mmax;i++)
{
tmp(i)=pow(v(i),x);
}
return tmp;
}
_______________________________________________
Users mailing list
Us...@admb-project.org
http://lists.admb-project.org/mailman/listinfo/users
dave fournier
Thanks Ben you found a fatal bug. Hard to believe that function had
never been used before. In the file f1b2vc1.cpp
the line
ADUNCONST(df1b2vector,x)
should have been
ADUNCONST(dvector,x)
df1b2vector operator + (const dvector& _x,const df1b2vector& _y)
{
ADUNCONST(df1b2vector,x)
ADUNCONST(df1b2vector,y)
check_shape(x,y,"operator +");
int mmin=x.indexmin();
int mmax=x.indexmax();
df1b2vector tmp;
tmp.noallocate(mmin,mmax);
for (int i=mmin;i<=mmax;i++)
{
tmp(i)=x(i)+y(i);
}
return tmp;
}
Sorry for that. The reason for these stupid errors is that every once
and a while the compiler police for the various c++ compilers make
changes that break huge amount of old code so one has to fix it all at
one go and these sorts of errors creep in.
Maybe someone can comment on how this can get fixed quickly.
as for just hoping that prob will remain between 0 and 1 -- that is a
mugs game. It will break at the most inopportune time. I modifed your
code to put in the bounds. also you need to include the prior for the u's
Finally WRT yuo commentabout initial values for the u's being necessary.
that is inaccurate. default values are provided for all parameters.
However if you want to supply a *.pin ifle then it must include values
for all parameters. I have already commented on the easiest way to
produce such a *.pin file.
Herre is your code as I have modified it.
//cout << "cvec calc\n";
cvec = sigma_c*(Z*u)+c;
a1vec = elem_prod(pow(size/d,g),exp(-size/d));
//cout << "prob calc\n";
prob = 1/(1/(elem_prod(cvec,a1vec))+h*initial);
dvariable ppen=0.0;
dvar_vector pprob=posfun(prob,0.01,ppen);
f+=ppen;
ppen=0.0;
dvar_vector pprob2=1.0-posfun(1.0-pprob,0.01,ppen);
f+=ppen;
f+=0.5*norm2(u); // prior
f -= sum( log_comb(initial,killed)+
elem_prod(killed,log(pprob2))+
elem_prod(initial-killed,log(1-pprob2)));
--
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com
Bolker,Benjamin Michael
Thanks very much Dave ...
Next question (sorry!!!)
Is there, somewhere, a start-up document about how to build
ADMB binaries from source? I say this because I feel like I am about
to enter dependency hell in trying to make/link the new version
of f1b2vc1.cpp -- e.g. when I try g++ -c f1b2vc1.cpp I get
warnings about .h files missing; when I copy the .h files over
from the ADMB source files, I then have to change the <file.h>
format to "file.h" so my system can find them ... etc. ...
It feels like the easiest (??) thing to do in order to
incorporate the bug fix that you sent would be to
replace the relevant file in the SVN copy of the source from
the google code site and rebuild the binaries, but that seems to contain "only" source code,
and not any information about the build environment (I guess I
was expecting to see Makefiles in there ...) Does anyone have
suggestions? I am running on Ubuntu Linux with standard tools (gcc/g++ 4.3.2 GNU Make 3.81, etc.)
Just doing something naive like cd admb-project-read-only/df1b2-separable; g++ -c *.cpp
doesn't work ...
I have poked around on
http://admb-project.org/community/public-domain
but not found anything I can understand yet ...
One of these days I will stop bothering everyone ...
cheers
Ben Bolker
________________________________________
From: users-...@admb-project.org [users-...@admb-project.org] On Behalf Of dave fournier [ot...@otter-rsch.com]
Sent: Thursday, August 06, 2009 11:29 AM
To: us...@admb-project.org
Subject: [ADMB Users] bug in file f1b2vc1
John Sibert
Below is a simple Makefile that works for the ADMB simple example with
Ubuntu (64 bit). You would have to modify it to include additional
source files. (But why don't you simply include the f1b2vc1.cpp source
code in the GLOBALS_SECTION of your tpl?).
Cheers,
John
(Remember the leading spaces need to be replaced by a tab.)
.SUFFIXES: .tpl
.PHONY: CLEAN RULES
CC=$(ADMB_HOME)/bin/mygcco
LL=$(CC)
simple: simple.cpp simple.o
$(LL) simple
.cpp.o:
$(CC) $*
.tpl.cpp:
$(ADMB_HOME)/bin/tpl2cpp $*
CLEAN:
rm -fv *.o
rm -fv *.htp
rm -fv simple.cpp
rm -fv simple
RULES:
@echo "ADMB_HOME = "$(ADMB_HOME)
@echo "CC = "$(CC)
@echo "LL = "$(LL)
--
Visit the ADMB project http://admb-project.org/
dave fournier
Just as before with the pow function you simply include the file
f1b2vc1.cpp into the GLOBALS_SECTION. You will have to add 2 spaces at
the beginning of each line because tpl2rem wants that. That's it!
Ben Bolker
df1b2s library, one gets lots of "previously defined" errors. If
I delete the -ldf1b2s flags from the compile statement:
g++ -s -L/usr/local/src/admb/lib mccoypred5.o -ladmod -ladt -lado -o
mccoypred5
then I get lots of stuff missing because there are many functions
in libdf1bs.a other than the new functions in f1b2vc1.cpp ...
Still trying to figure out the build process ...
thanks,
Ben
dave fournier wrote:
> Sorry I thought it was clear what to do.
> Just as before with the pow function you simply include the file
> f1b2vc1.cpp into the GLOBALS_SECTION. You will have to add 2 spaces at
> the beginning of each line because tpl2rem wants that. That's it!
>
--
Ben Bolker
Associate professor, Biology Dep't, Univ. of Florida
bol...@ufl.edu / www.zoology.ufl.edu/bolker
GPG key: www.zoology.ufl.edu/bolker/benbolker-publickey.asc