Strange State Assignment estimates and EM(20) not working

[simone77]

Hi Rémi,

Thanks. Regarding:

However, I think that you should rewrite the model for survival.

Let consider two shortcuts:

Juveniles = [g(1:3).a(1)]
Adults = [g(1:3).a(2)+g(4:6)]

A simple model with two "real classes of age" is

Juveniles+Adults

A model with the same effects used for your general model is

Juveniles.f+Adults.f.t

Yes! I was forgetting juveniles group become adults after first capture. Haste makes waste...

[simone77]

Just to avoid confusion for someone else reading this thread: I believe there was an error in the shortcut for adults. Instead of:
Adults = [g(1:3).a(2)+g(4:6)]
It would be:
Adults = [g(1:3).a(2)&g(4:6)]

[jbleu]

Hi,

I am using E-surge 1.8.5
I am working on a dataset of females and I am interested in the influence of reproduction on survival. Thus I have 3 states (reproduced, did not reproduce, dead) and 4 events (not seen, seen and she reproduced, seen and did not reproduce, seen but reproductive state is undetermined). I have 21 years of data.
I created the encounter histories with 0,1,2 and 3 corresponding to each event.

In GEPAT I have:
%%%% Initial state %%%%%%
1
1 2 IS
p *
%%%% Transition %%%%%%
2
3 3 S (survival probability)
y - *
- y *
- - *
3 3 R (transition between reproductive states)
y * -
y * -
- - *
%%%% Event %%%%%%
2
3 3 C (probability to capture/observe the females)
* b -
* - b
* - -
3 4 SA (state assignment)
* - - -
- b - *
- - b *

and in Gemaco I implemented a simple model:
IS=i
S=from
R=from
C=t
SA=from

and I have the same problem as Simone:

1. EM(20)+Quasi-newton not working (see above first post)
2. SA (State Assignment) estimates for one gender are nonsensical: are all "1s".

In my case, SA estimates the probability to assign the reproductive state to 1 for non reproductive females and to 0.67 for reproductive females... This does not make a lot of sense that the females in the "unkown" state are all reproductive females...
(If I model a constant SA, I obtain 0.72).

I did not understand how Simone solved this problem of strange state assignment...
Any suggestions?

Thank you for your help
Josefa

[simone77]

Hi Josefa,

Based just on my experience I can tell you the EM(20) Quasi-Newton procedure doesn't work always. I have been trying it on several different analyses (exploratory or not) and sometimes some error warning appeared and it didn't work.
Remarkably, as Rémi told in this topic, it is a procedure good to deal with boundary estimates but I have no idea about the reason why it sometimes works and others doesn't.

About your case of study, you have two steps for events modeling (i) the p to be seen and (ii) the p the reproductive status is assessed (you assume there is no error when breeding status is assessed). I am not sure I understood your interpretation of the results: did you get estimates of SA from (2) to (2) = 0.67 and from (3) to (3) = 1?
If it were so, according to your GEPAT structure, the reproductive status of a breeding female is assessed 67% of times whereas that of a non breeding female is always assessed (p=1).
In other words, a possible explanation might be that non-reproductive females are very easy to be assessed as non-reproductive females whereas the reproductive status of breeding females may be assessed or not. Does it make sense?

I apologize if I misunderstood something.

[jbleu]

Thank you for your quick reply.

I was afraid that if the EM(20) Quasi-Newton doesn't work it meant that my model was not correctly defined in gepat or gemaco. If it is not supposed to work all the time, then it is not really a problem for me.

About your case of study, you have two steps for events modeling (i) the p to be seen and (ii) the p the reproductive status is assessed (you assume there is no error when breeding status is assessed).

yes exactly.

did you get estimates of SA from (2) to (2) = 0.67 and from (3) to (3) = 1?

yes this is correct.

If it were so, according to your GEPAT structure, the reproductive status of a breeding female is assessed 67% of times whereas that of a non breeding female is always assessed (p=1).
In other words, a possible explanation might be that non-reproductive females are very easy to be assessed as non-reproductive females whereas the reproductive status of breeding females may be assessed or not. Does it make sense?

I have the same interpretation, but it does not make biological sense.
Females are classified "reproduced" if they are observed with a young. however, females are not always with their young and thus if the female is not observed very long we don't know if she didn't reproduce or if she is just away from her kid. Thus the female is classified as "unknown status". The females are also classified "unknown" if they are observed early in the season and thus we don't know if she reproduced or not.
Thus I don't expect that reproducing or non reproducing females are always assessed. A p=1 is not biologically meaningful...

Also, I tried to run the model with multiple random initial values, but I have the same issue...

Thanks for your help

[simone77]

Hi Josefa,

Definitively running models with multiple random (MR) initial values is always a good idea because it is a good way to deal with local minima. Personally, in the exploratory phase I use to run the same model with different numbers of MR initial values just to have an idea on how local minima are present or not and what is the minimum threshold number of MR needed to "get rid" of them. Another useful trick is to retrieve the model and check the "from last model" option and run it again once it is retrieved (you must not to set anything in the IVFV option in this case), sometimes you will get a lower deviance from this model. All this stuff is useful to have an idea about the quality of your data multi-event analysis and sometimes they allow you to get an estimate that otherwise was too close to the boundary and estimated as "1".

That said, I am a bit confused about your analysis: when a female is classified as non-reproductive?

[jbleu]

Hi,
We have observations of the females during spring and summer, however the reproductive state can only be assessed during summer. So if a female is seen only during spring we don't know the reproductive state. A female is non-reproductive only when she is seen during summer without a kid (females produce only one kid) AND that she is doing an activity where normally the females do not leave their kids (for example it is frequent that females forage alone but when they are resting we expect to see the kid nearby for a reproducing female). Also if the female is observed during a long time (several hours) without a kid she is non-reproductive.
In the dataset we have: number of reproductive females > number of unknown females > number of non-reproductive females.
Can it be a problem for state assignment that the number of undetermined females is higher than the number of non-reproducing females ?

Josefa

[simone77]

Hi Josefa,

I have to highlight that I am just an user (like you) and therefore I can be missing the point, so take it into consideration...

A priori, if the number of known non-breeding females is not excessively low, I wouldn't expect that its relative number with respect to the known breeding females and the breeding unknowns are beyond these results. Anyway, I would have a look at the raw numbers to have an idea about how much sparse is the information you are trying to modeling. Have you run before U-Care on your data to check the GOF of the CJS (often if the number of observations with known states are low, the unistate CJS model is tested)?

One possibility would be create a fake data set with a higher number of "2's" (event for non-breeding females) and repeat the same identical analysis on it and see what happens: do you still have "1" for SA estimate of non-breeding females? if so, it might suggest something is wrong with the structure of your analysis more than with the quality of your data.

[jbleu]

Hi,

I had not run U-care on my data. Thank you for the advice !
Now, I have done GOF in U-CARE on unistate encounter histories on the whole data set. TEST3.SR + TEST3.SM + TEST2.CL was no significant, thus it suggests to use as global model the Phi(t) p(t*m) (with trap dependency). Also if I test the jmv model with the multistate model I have: TEST WBWA + TEST 3G.SR + TEST 3G.SM + TEST M.LTEC not significant.

Thus I created an umbrella model, with a gemaco like this :
IS = i
S= from.t
transitions between state = from.t
capture = firste + nexte.from.t
state assignment = from

also ivfv of first "event parameter" is set to 1 (first capture probability).

The results for SA are from (2) to (2) = 0,99999415 and from (3) to (3) = 0,342042085
I still have a probability of SA close to 1 but it is the opposite than in the previous model (in the previous model from (3) to (3) was equal to 1).
Thus it seems that depending on model definition one of the two probability is 1, but not always the same...

Also during the iteration, this sentence appeared: "10 first histories incompatible with the model". I don't know what is wrong with these histories. Do you know what I should check ? Is it really a problem: it did not prevent the iteration to converge, and (i think that) it did not appear when I ran the model a second time with different initial values...

thanks for your help !

[simone77]

Hi,

I had not run U-care on my data. Thank you for the advice !
Now, I have done GOF in U-CARE on unistate encounter histories on the whole data set. TEST3.SR + TEST3.SM + TEST2.CL was no significant, thus it suggests to use as global model the Phi(t) p(t*m) (with trap dependency). Also if I test the jmv model with the multistate model I have: TEST WBWA + TEST 3G.SR + TEST 3G.SM + TEST M.LTEC not significant.

Trap dependence (trap happiness or shyness) can be dealt in several ways: (i) you can modify the data set and use U-Care to do that (Pradel R (1993) Flexibility in survival analysis from recapture data: handling trap-dependence. In: Lebreton J-D, North PM (eds) Study of bird population. Birkhauser Verlag, Basel, pp 29–37), (ii) or, as proposed recently by Pradel and Sanz-aguilar, modify the GEPAT structure of your analysis (Pradel & Sanz-Aguilar 2012), or (iii) simply give a different proability of resighting for the first interval after the first capture (the first and latter are discussed here).

Thus I created an umbrella model, with a gemaco like this :
IS = i
S= from.t
transitions between state = from.t
capture = firste + nexte.from.t
state assignment = from

also ivfv of first "event parameter" is set to 1 (first capture probability).

Have you any biological reason to believe that the probability one first-captured individual were a breeder is constant across the capture sessions? if you have not, I believe it would more appropiate to let IS to depend on unspecific time variation. Also, I don't believe that you are modeling trap-dependence with that syntaxis on Capture (see the above comments).

The results for SA are from (2) to (2) = 0,99999415 and from (3) to (3) = 0,342042085
I still have a probability of SA close to 1 but it is the opposite than in the previous model (in the previous model from (3) to (3) was equal to 1).
Thus it seems that depending on model definition one of the two probability is 1, but not always the same...

It is not easy to me to find an explication about this from the distance, it seems like something is going wrong with the way the data set is built or perhaps in the GEPAT....really don't know where the problem can arise.

Also during the iteration, this sentence appeared: "10 first histories incompatible with the model". I don't know what is wrong with these histories. Do you know what I should check ? Is it really a problem: it did not prevent the iteration to converge, and (i think that) it did not appear when I ran the model a second time with different initial values...

thanks for your help !

That's not necessarily a problem, not if the model converges. Remi Choquet answered this same question (from me!) here.

Hope this may help.
Good luck!