Negative AIC

[Diego.Pavon]

Hello,

I am running models (type: psi() eps() p()) with food abundance as a covariate (I entered directly in the DM). My problem is that in the ranking, my first model (the best?) has AIC = -1; DeltaAIC= -54 and weight = 1. This model is: psi(prey)eps(prey)p(seasonal effect). My second model in that ranking (psi(prey)eps(.)p(seasonal effect+prey)) has AIC = 53, deltaAIC = 0 and weight = 0. The 3rd model has already a deltaAIC = 1200...

What does it mean? is the first model with negative AIC the best? Does it mean that something when calculating AIC went wrong?

Thank you

Diego

[jhines]

Diego,

I think the AIC=-1 means something is wrong. Can you email me (jhines@usgs.gov) the latest presence_backup.zip file in your project folder?

Jim

[jhines]

Diego,

The model with AIC=-1 did not converge. If you look at the output from that model, you'll see the value for AIC is -1.#IND00, which means the value could not be computed due to an arithmetic problem (divide by zero). This model should be deleted.

A problem with the psi(t),eps(t) parameterization is that since gamma is not estimated, but computed from psi and eps, it's possible for the computed gamma to be < zero or > 1. When this happens, some of the computations are impossible and the likelihood function becomes that value you saw in the AIC (-1.#IND00). A possible solution is to give PRESENCE better initial values to start the iterative process so that doesn't happen. One of your models (psi()eps(voles)p(voles) is nearly identical to the problem model (psi(voles)eps(voles)p(voles)), but it did converge. I suggest taking the final beta estimates from that model and using them as starting values for the problem model. Since the other model only estimates one beta for psi, you'll have to insert a zero for the 2nd psi beta. To give starting values for a model, check the box labeled 'specify starting values' in the run model dialog box. When the box appears, enter the final beta estimates from the previous model into the box (one entry per line), then insert a line with a '0' just after the 1st line.

I think a bigger concern here is that you're putting structure on two parameters that depend on each other. Seasonal occupancy in time t+1 depends on extinction in time t.

psi(t+1) = psi(t)*(1-eps(t)) + (1-psi(t))*gam(t)

You can think of seasonal occupancy as being a function of voles, or extinction being a function of voles, but not both. Actually, if extinction is a function of voles, occupancy will also be dependent on voles, since occupancy depends on extinction. I suggest only running models with the covariate on occupancy and allowing extinction to be time-specific, or models with occupancy time-specific and extinction as a function of voles.

One last thing, it looks like the extinction values in the output are really persistence (1 - extinction). This is a bug in the program and will be corrected soon.

Cheers,

Jim

[Diego.Pavon]

Hi Jim,

Well, it does not sound good. Let see if I understood your suggestion:

If I want to enter initial values, I should enter one value for each psi, eps and p, am I right? That is, if I have 22 years (22 ocassions), and I want to enter initial values for the model psi(t)eps(voles)p(.) I need to enter: 22 values (rows) for psi, 2 values for eps (one for intercept and one for the covariate,i.e. voles) and 1 value for p (the intercept). Is like this how it should be?

Another thing that comes to my mind is whether or not the estimates that I will obtain from these models depend on the initial values. I mean, Would the same model give me different estimates of psi,eps and p, depending on what are the initial values?
In that sense, how to choose the initial values? You suggested that the best way is to pick them from a similar model, but if there are several "similar" models, should one try them all?

I have found that some other models have not converged. Maybe that is why it took me more than 8 hours to run some models...

Thank you,

Diego

[jhines]

Diego,

Yes, you need to enter 22 initial values for psi, as there are 22 columns in the first design matrix. Two initial values for eps (1 intercept, 1 covariate effect), and one initial value for p. You don't need to be too precise with the initial values, and I've found that the most important initial values are the ones related to parameters which depend on covariates. So, if the beta estimate for occupancy from a model with psi constant is .25, I'd enter .25 22 times. If you think extinction probability is relative low, I'd enter -2 for the intercept, and if you think voles might cause your species to decline, I'd enter 1.0 for the slope or effect parameter (thinking higher number of voles leads to higher extinction). Finally, you could enter the last beta estimate from another model which has p constant for the last initial value.
I wouldn't worry too much about trying lots of initial values. Usually, the optimization routine in PRESENCE is pretty good at finding the global maximum likelihood value, regardless of initial values. The cases where it has trouble are when data are sparse and/or the model is overparameterized (ie., not enough data or trying to estimate too many parameters). Signs of this in the output are: Error or Warning messages, extremely high standard errors for the estimates, or Log-likelihood or AIC values which are extremely different from other models. If you don't see any of those signs, then the model most-likely converged on the global maximum-likelihood value and I wouldn't worry about trying other initial values. (In the model you initially asked about, all of these signs were there.)
A couple of exceptions to this recommendation are: 1) Warning about non convergence and 2) high standard errors of estimates which are <-10 or >10. Due to computer rounding errors, the optimization routine will sometimes report that the model has not fully converged and is only accurate to X.XX significant digits. If this number (X.XX) is greater than 2, I’d ignore the warning as the estimates will be accurate to at least 4 digits. Estimates for beta parameters which correspond to real parameters which are near a boundary (zero or one), will usually have high standard errors. This is expected as the standard error of a real parameter which is zero or one is undefined. So, those standard errors can be ignored.
If you see signs of problems with convergence, then I recommend trying other initial values. If the log-likelihood value from the output with different initial values is higher than the output in question, then you can delete the original output and use the one with initial values. (The last column in the model summary screen in PRESENCE gives -2*LogLike, so you’d look for a model with a lower value in this column.) I think the models you’re trying are not so complicated that you need to try a lot of initial values. Also, the log-likelihood values are pretty consistent among models.
Jim

[Diego.Pavon]

Hi Jim,

Thank you for the explanation! It is great. I have some doubts, though, about when to keep a model as a good model: If I understood correctly, if the warning message is, e.g.:
**** Numerical convergence was not reached.
Parameter estimates converged to approximately 4.52 significant digits.

I should keep the model (the number is > 2).

However, the estimates (betas) are not so good looking (for instance):
estimate std.error
A1 ccupancy psi1 0.671743 (6687834.736011)
A2 ccupancy psi1 0.079996 (10352644.074366)
C1 :local extinction eps1 1.346402 (1.553351)
C2 :local extinction eps2 2.356490 (1.327408)
C3 :local extinction eps3 24.510901 (0.000011)
C8 :local extinction eps8 36.048295 (0.000000)
C14 :local extinction eps14 24.716508 (-1.#IND00)

Should I still keep this model as a correct one? In this case it is not so important because it is far from the best model (>10 AIC units) but I am wondering whether this model (if converged) would be up in the ranking. I have gotten many models with this warning message with significant digits >2, but the errors and estimates are a bit weird.
Shall I then keep trying until I get no warning message at all?

Another thing: You suggested not to use the same variable for psi and extinction. I agree. But the thing is that, in case of voles, I constrain Psi as a function on voles in preceding autumn, that is, PSI in season t is function of vole abundance in t - 1 (the same applies to p). Instead, epsilon is function of the voles in that year (e.g. Epsilon3 which is the probability to get extinct between season 3 to season 4 is function of vole abundance in season 3). Therefore, there are "different" variables. Do you think that I should not have this model: PSI(voles preceding autumn)EPS(voles autumn)? Do you still think that I should have only one of them as a function of voles in the same model?

Thank you,
Diego

[jhines]

Diego,

If the log-likelihood value from that model is not radically different from the other models, then I would not bother trying other starting values. I suspect that it is overparameterized and has no chance of being a top model. It doesn't matter if you delete it or keep it as it will have a model weight of zero.

If this model was one of your top models, then I would try to resolve the problems with the variances. I usually start by looking at the estimates with very large variances (or variance of zero) to see if there is a reason why they cannot be estimated. If those estimates are on a boundary (logit transformed value of zero or one), then I would fix the real parameter associated with it to zero or one, and delete that column from the design matrix.

With overparameterized models, different starting values will not help in getting reasonable estimates and variances. The model has to be simplified in some way until the model is not overparameterized.

Again, I would not bother with this model since it's not going to be near the top model.

Jim