May be a problem with Huggins closed on latest version?

[DuaneD]

I repeated an analysis of a removal model (3 occasions) where I used length as an icovariate. I fixed recapture probably (c) = 0. In the design matrix I have an intercept and slope parameter for length.

In the latest version of MARK (downloaded and installed today 2/11/15), the program estimates a beta parameter for recapture probability (even though fixed) and I get a message that there is an estimability problem. However, the abundance estimate matches the previous analysis.

Below I paste the model output from 2005 and then today's analysis. Let me know if you'd like the input or MARK project files.

Duane


Program MARK - Survival Rate Estimation with Capture-Recapture Data
Version 4.2(Win32) Apr 2005 14-Oct-2005 12:18:30 Page 001
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -



INPUT --- proc title Huggins Removal Model for Brown Trout, Spruce
INPUT --- Creek, Pennsylvania;

Time in seconds for last procedure was 0.00


INPUT --- proc chmatrix occasions=3 groups=1 etype=Huggins icovar=2
INPUT --- Standardize ICMeans NoHist hist=300;

INPUT --- glabel(1)=Trout;

INPUT --- time interval 1 1 1;

INPUT --- icovariates age length;

Number of unique encounter histories read was 49.

Number of individual covariates read was 2.
Time interval lengths are all equal to 1.

Data type is Huggins closed captures.

Time in seconds for last procedure was 0.03

 Program MARK - Survival Rate Estimation with Capture-Recapture Data
Version 4.2(Win32) Apr 2005 14-Oct-2005 12:18:30 Page 002
Huggins Removal Model for Brown Trout, Spruce Creek, Pennsylvania
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -



INPUT --- proc estimate link=Logit varest=2ndPart ;


INPUT --- model={p(length)};

INPUT --- fixed=2;
INPUT --- parm(4)=0;
INPUT --- parm(5)=0;

INPUT --- group=1 p rows=1 cols=3 Square;
INPUT --- 1 2 3;

INPUT --- group=1 c rows=1 cols=2 Square;
INPUT --- 4 5;

INPUT --- design matrix constraints=5 covariates=3;
INPUT --- 1 length 0;
INPUT --- 1 length 0;
INPUT --- 1 length 0;
INPUT --- 0 0 1;
INPUT --- 0 0 1;

INPUT --- blabel(1)=;

INPUT --- blabel(2)=;

INPUT --- blabel(3)=;

INPUT --- rlabel(1)=p;

INPUT --- rlabel(2)=p;

INPUT --- rlabel(3)=p;

INPUT --- rlabel(4)=c;

INPUT --- rlabel(5)=c;

Link Function Used is LOGIT

Variance Estimation Procedure Used is 2ndPart

 Program MARK - Survival Rate Estimation with Capture-Recapture Data
Version 4.2(Win32) Apr 2005 14-Oct-2005 12:18:30 Page 003
Huggins Removal Model for Brown Trout, Spruce Creek, Pennsylvania
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -


M(t+1):
747


-2logL(saturated) = 0.0000000
Effective Sample Size = 2241

Number of function evaluations was 12 for 2 parameters.
Time for numerical optimization was 0.02 seconds.
-2logL {p(length)} = 737.17546
Penalty {p(length)} = 0.0000000
Gradient {p(length)}:
-0.4358380E-05 0.1160941E-04
S Vector {p(length)}:
93.49314 30.64193
Time to compute number of parameters was 0.01 seconds.
Threshold = 0.6000000E-07 Condition index = 0.3277453
Conditioned S Vector {p(length)}:
1.000000 0.3277453
Number of Estimated Parameters {p(length)} = 2
DEVIANCE {p(length)} = 737.17546
DEVIANCE Degrees of Freedom {p(length)} = 745
c-hat {p(length)} = 0.9894973
AIC {p(length)} = 741.17546
AICc {p(length)} = 741.18082
Pearson Chisquare {p(length)} = 14908.624


LOGIT Link Function Parameters of {p(length)}
95% Confidence Interval
Parameter Beta Standard Error Lower Upper
------------------------- -------------- -------------- -------------- --------------
1: 1.6084657 0.1042736 1.4040894 1.8128421
2: 0.9585298 0.1801610 0.6054142 1.3116455
3: 0.0000000 0.0000000 0.0000000 0.0000000


Real Function Parameters of {p(length)}

Following estimates based on standardized individual covariate values:
Variable Value Mean SD
--------- ------------- ------------- -------------
AGE 0.0 0.6530612 0.4809288
LENGTH 0.0 234.30612 57.666500
95% Confidence Interval
Parameter Estimate Standard Error Lower Upper
------------------------- -------------- -------------- -------------- --------------
1:p 0.8331983 0.0144918 0.8028320 0.8597050
2:p 0.8331983 0.0144918 0.8028320 0.8597050
3:p 0.8331983 0.0144918 0.8028320 0.8597050
4:c 0.0000000 0.0000000 0.0000000 0.0000000 Fixed
5:c 0.0000000 0.0000000 0.0000000 0.0000000 Fixed

 Program MARK - Survival Rate Estimation with Capture-Recapture Data
Version 4.2(Win32) Apr 2005 14-Oct-2005 12:18:30 Page 004
Huggins Removal Model for Brown Trout, Spruce Creek, Pennsylvania
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -



Estimates of Derived Parameters
Population Estimates of {p(length)}
95% Confidence Interval
Group N-hat Standard Error Lower Upper
--------- -------------- -------------- -------------- --------------
1 754.06921 3.5510479 749.78995 764.91203

Time in seconds for last procedure was 0.03


INPUT --- proc stop;

Time in minutes for this job was 0.00


E X E C U T I O N S U C C E S S F U L


Here's the current run:


Program MARK - Survival Rate Estimation with Capture-Recapture Data
gfortran(Win64) Vers. 8.0 Jun 2014 11-Feb-2015 13:39:08 Page 001
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

This version was compiled by GCC version 4.8.1 using the options
-m64 -mtune=generic -march=x86-64 -mthreads -O2 -fimplicit-none -fbounds-check -funroll-loops -ftree-vectorize -fopenmp.
This problem will use 8 of 8 possible threads.


INPUT --- proc title f;

CPU Time in seconds for last procedure was 0.00


INPUT --- proc chmatrix occasions=3 groups=1 etype=Huggins icovar=2
INPUT --- ICMeans NoHist hist=300;

INPUT --- glabel(1)=Group 1;

INPUT --- time interval 1 1 1;

INPUT --- icovariates age len;

Number of unique encounter histories read was 49.

Number of individual covariates read was 2.
Time interval lengths are all equal to 1.

Data type number is 12
Data type is Huggins' p and c

CPU Time in seconds for last procedure was 0.00

 Program MARK - Survival Rate Estimation with Capture-Recapture Data
gfortran(Win64) Vers. 8.0 Jun 2014 11-Feb-2015 13:39:08 Page 002
f
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -



INPUT --- proc estimate link=Logit varest=2ndPart ;


INPUT --- model={p(len)};

INPUT --- fixed=2;
INPUT --- parm(4)=0;
INPUT --- parm(5)=0;

INPUT --- group=1 p rows=1 cols=3 Square Time=1;

INPUT --- group=1 c rows=1 cols=2 Square Time=4;

INPUT --- design matrix constraints=5 covariates=3;
INPUT --- 1 len 0;
INPUT --- 1 len 0;
INPUT --- 1 len 0;
INPUT --- 0 0 1;
INPUT --- 0 0 1;

INPUT --- blabel(1)=;

INPUT --- blabel(2)=;

INPUT --- blabel(3)=;

INPUT --- rlabel(1)=p;

INPUT --- rlabel(2)=p;

INPUT --- rlabel(3)=p;

INPUT --- rlabel(4)=c;

INPUT --- rlabel(5)=c;

Link Function Used is LOGIT

Variance Estimation Procedure Used is 2ndPart

M(t+1):
747



 Program MARK - Survival Rate Estimation with Capture-Recapture Data
gfortran(Win64) Vers. 8.0 Jun 2014 11-Feb-2015 13:39:08 Page 003
f
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

-2logL(saturated) = -0.0000000
Effective Sample Size = 2241

Number of function evaluations was 10 for 3 parameters.
Time for numerical optimization was 0.06 seconds.
-2logL {p(len)} = 737.17546
Penalty {p(len)} = -0.0000000
Gradient {p(len)}:
0.000000 0.1631725E-05-0.1033546E-04
S Vector {p(len)}:
131.6585 0.5614433 0.1017915E-04
Time to compute number of parameters was 0.01 seconds.
Threshold = 0.8000000E-07 Condition index = 0.7731479E-07 New Threshold = 0.1763080E-04
New Guessimate of Estimated Parameters {p(len)} = 2
Conditioned S Vector {p(len)}:
1.000000 0.4264391E-02 0.7731479E-07
Number of Estimated Parameters {p(len)} = 2
DEVIANCE {p(len)} = 737.17546
DEVIANCE Degrees of Freedom {p(len)} = 745
c-hat {p(len)} = 0.9894973
AIC {p(len)} = 741.17546
AICc {p(len)} = 741.18082
BIC {p(len)} = 752.60481
Pearson Chisquare {p(len)} = 1761.3121


LOGIT Link Function Parameters of {p(len)}
95% Confidence Interval
Parameter Beta Standard Error Lower Upper
------------------------- -------------- -------------- -------------- --------------
1: -2.2861593 0.7285402 -3.7140980 -0.8582205
2: 0.0166220 0.0031242 0.0104985 0.0227454
3: -0.0999689 0.0000000 -0.0999689 -0.0999689


Real Function Parameters of {p(len)}

Following estimates based on unstandardized individual covariate values:
Variable Value
--------- -------------
AGE 0.8206158
LEN 240.55288
95% Confidence Interval
Parameter Estimate Standard Error Lower Upper
-------------------------- -------------- -------------- -------------- --------------
1:p 0.8471342 0.0139948 0.8176464 0.8725967
2:p 0.8471342 0.0139948 0.8176464 0.8725967
3:p 0.8471342 0.0139948 0.8176464 0.8725967
4:c 0.0000000 0.0000000 0.0000000 0.0000000 Fixed
5:c 0.0000000 0.0000000 0.0000000 0.0000000 Fixed

 Program MARK - Survival Rate Estimation with Capture-Recapture Data
gfortran(Win64) Vers. 8.0 Jun 2014 11-Feb-2015 13:39:08 Page 004
f
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -



Estimates of Derived Parameters
Population Estimates of {p(len)}
95% Confidence Interval
Group N-hat Standard Error Lower Upper
--------- -------------- -------------- -------------- --------------
1 754.06921 3.5513623 749.78975 764.91335

Attempted ordering of parameters by estimatibility:
1 2 3
Beta number 3 is a singular value.

CPU Time in seconds for last procedure was 0.09


INPUT --- proc stop;

CPU Time in minutes for this job was 0.00

Time Start = 13:39:08.871 Time End = 13:39:08.918

Wall Clock Time in minutes for this job was 0.00


E X E C U T I O N S U C C E S S F U L

[cooch]

Hi Duane --

Quick answer to some of your query. I'm already a bit late for something.

I repeated an analysis of a removal model (3 occasions) where I used length as an icovariate. I fixed recapture probably (c) = 0. In the design matrix I have an intercept and slope parameter for length.

In the latest version of MARK (downloaded and installed today 2/11/15), the program estimates a beta parameter for recapture probability (even though fixed) and I get a message that there is an estimability problem. However, the abundance estimate matches the previous analysis.

For a fixed parameter, MARK isn't actually estimating anything -- there is a beta listed in the output for the fixed parameter, but it is a nonsense value, usually the transform of the MARK default starting parameter 'guess' (typically 0.45 for parameters on the [0,1] interval).

For example, take the following removal data set:

Code: Select all

100000 52;
010000 19;
001000 14;
000100 7;
000010 5;
000001 2;   


If I fit a Huggins model p(t)c(0) (i.e., where c is constant, fixed to 0), using the logit link, then my beta estimates are:

Code: Select all

                  LOGIT Link Function Parameters of {phi(t)p(0)} - logit}
                                                              95% Confidence Interval
 Parameter                    Beta         Standard Error      Lower           Upper
 ----------------------  --------------  --------------  --------------  --------------
    1:p                     0.0538529       40.140183       -78.620907      78.728612
    2:p                     -0.4658409      65.332352       -128.51725      127.58557
    3:p                     -0.1504940      121.53586       -238.36079      238.05980
    4:p                     -0.2812736      213.27325       -418.29684      417.73430
    5:p                     0.1569540       462.78948       -906.91045      907.22436
    6:p                     -0.1282681      869.86859       -1705.0707      1704.8142
    7:c                     -0.1000032      0.0000000       -0.1000032      -0.1000032       


So, yes, there is something reported for parameter c, but I wouldn't go so far as to say it is estimated (a 0 SE and UCI/LCI being equal is pretty diagnostic).

And -0.1 is roughly logit(0.45). Close enough for government work. Whoops -- you're a Fed.

[cooch]

Example of the same thing...multi-state, two states (A, B), 5 occasions. No marked individuals in state A on occasion 1. Since there are no individuals in state A, then the transition from A -> B over first interval can't be estimated. So, if * *don't* fix psi to be a logical 0, MARK will happily report the real estimate for that parameter as 0.45, which of course is the MARK default starting value. The reported SE of the 'estimate' is 0. Meaning, in this case, nothing gets estimated.

If you look at the beta's, you'll see -0.1000. Now what is potentially somewhat confusing is that what MARK reports for the beta for these logical 0's is more or less the same value (-0.1000) regardless of the link function you use (although they do differ somewhat in the 5th decimal place and beyond). Back-transform from -0.1 on the sin scale is ~0.45, whereas on the logit scale, it 0.475, which is what MARK reports.

I'll leave it to Gary to explain that nuance -- my *guess* is that Gary took the default link (sin), figured out the sin transformed value for a real starting value of 0.45 (which is -0.10016742....), and then simply used that numerical value as the starting point for all of the different link functions. Again, just a guess.

[cooch]

I'll leave it to Gary to explain that nuance -- my *guess* is that Gary took the default link (sin), figured out the sin transformed value for a real starting value of 0.45 (which is -0.10016742....), and then simply used that numerical value as the starting point for all of the different link functions. Again, just a guess.

Nope -- must be something else. I just tried the log link, and the estimated beta (-0.798413) is 0.45 when back-transformed.

So, for some reason which is probably unimportant, MARK seems to use 0.475 (on the real scale) as the starting value for models estimated using the logit link, and 0.45 for everything else.

[DuaneD]

Ok, I see now that when you fix parameters MARK no longer removes the parameter from the list to be estimated because you now get the statement in the model run

"Number of function evaluations was 10 for 3 parameters." even when 1 is fixed,

whereas way back in 2005 the same model run stated

"Number of function evaluations was 12 for 2 parameters."

I guess I need to ignore estimability problems when they refer to a fixed parameter (and assume there are no other issues).

Sorry for any inconvenience,

Duane

[cooch]

Ok, I see now that when you fix parameters MARK no longer removes the parameter from the list to be estimated because you now get the statement in the model run

"Number of function evaluations was 10 for 3 parameters." even when 1 is fixed,

whereas way back in 2005 the same model run stated

"Number of function evaluations was 12 for 2 parameters."

I guess I need to ignore estimability problems when they refer to a fixed parameter (and assume there are no other issues).

Still a fair question as to why Gary made the change (presumably, on purpose). Reporting the fixed parameters in the list of betas, even if the 'estimates' are clearly non-sensical, has the potential to be confusing. Unless there is an obvious benefit to reporting them, I'd be inclined to ask Gary if it might not be worth 'changing back'.

[gwhite]

First, I'm not sure why you got an estimation issue, because I don't normally see that when you fix the parameters. I don't see any warnings in the output, so I'm wondering about that.

There are at least 2 things you can do to avoid any issues with the beta estimates when fixing real parameters. First just zero out the DM columns that pertain to the fixed parameters. MARK checks for DM columns with all zeros, and just ignores them. So in your case, the DM would look like
INPUT --- design matrix constraints=5 covariates=3;
INPUT --- 1 length 0;
INPUT --- 1 length 0;
INPUT --- 1 length 0;
INPUT --- 0 0 0;
INPUT --- 0 0 0;

However, even better would be to just delete the whole column, because it is not needed:

INPUT --- design matrix constraints=5 covariates=2;
INPUT --- 1 length;
INPUT --- 1 length;
INPUT --- 1 length;
INPUT --- 0 0;
INPUT --- 0 0;

Gary