F.cjs.estim(capture, survival, histories, cap.init, sur.init, group,
nhat.v.meth = 1, c.hat = -1, df = NA,
intervals=rep(1,ncol(histories)-1), conf=0.95, control=mra.control())capture. The default value usually works.survival. The default value usually works.nhat.v.meth = 1 uses the variance estimator of
Taylor et al. 2002, Ursus, p. 188 which is the so-called Huggins variance estimator, and incorporatec.hat) to use
during estimation. If input value of c.hat is df == NA, the number of parameters is estimated from the rank of the
matrix of 2nd derivatives or Hessian, depending on cov.meth parameter.
If <ncol(histories)-1 (i.e.,
number of capture occasions minus 1) specifying relative time intervals between occasions.
For example, if capture occasions occurred in 1999, 2000, 2005, and 2007 print.cr to print
it nicely. Use names(fit), where the call was fit <- F.cr.estim(...),
to see names of all returned components. To see values of individual components,
issue commands like fit$s.hat, fit$se.s.hat, fit$n.hat, etc.
Components of the returned object are as follows:loglik. This is relative deviance, see help for
F.sat.lik.deviance + 2*(df). df is either the estimated number of independent
parameters (by default), or NX+NY, or a specified value, depending on the input value of DF parameter.deviance / vif) + 2(df)df*(df+1)) / (nan - df - 1)df*(df+1))/(nan - df - 1)chisq.vif / chisq.dfvif, based
on pooling rules.exit.code), and
(3) covariance matrix code (interprets cov.code).exit.code is in message.exit.code = 3, fn.evals equals the maximum set in mra.control. This, in combination
with the exit codes and execution time, can help detect non-convergence or bad behavior.exit.code, cov.code, and fn.evals, this could be used to detect ill-behaved
or marginally unstable problems, if you know what you are doing. Assuming maxfn is set high
in mra.control() (e.g., 1000),
if exit.code = 1 but the model takes a long time to execute relative to similarly sized problems,
it could indicate unstable or marginally ill-behaved models.n.hat = NS. No estimate for
first occasion.n.hat estimates. Computed using method
specified in nhat.v.meth.n.hat.conf percent on n.hat. Length NS.n.hat.conf percent on n.hat. Length NS.n.hatn.hathistories matrix.F.cjs.gof for a description of $\Psi_{ij}$.mra.log, is written to the current directory. This file contains
additional details, such as individual Test 2 and Test 3 components, in a semi-friendly
format. This file is overwritten each run.
Model Specification: Both the capture and survival model can be
specified as any combination of 2-d matrices (time and individual varying covariates),
1-d time varying vectors, 1-d individual
varying vectors, 1-d time varying factors, and 1-d individual varying factors.
Specification of time or individual varying effects uses the
tvar (for 'time varying') and ivar (for 'individual varying') functions.
These functions expand covariate vectors along the appropriate dimension to
be 2-d matrices suitable for fitting in the model. ivar expands
an individual varying vector to all occasions. tvar expands a time
varying covariate to all individuals. To do the expansion, both tvar
and ivar need to know the size of the 'other' dimension. Thus, tvar(x,100)
specifies a 2-d matrix with size 100 by length(x).
ivar(x,100)
specifies a 2-d matrix with size length(x) by 100.
For convenience, the 'other' dimension of time or individual varying covariates
can be specified as an attribute of the vector. Assuming x
is a NS vector and the 'nan' attribute of x has been set as
attr(x,"nan") <- NAN, tvar(x,NAN) and
tvar(x) are equivalent. Same, but vise-versa, for individual varying
covariates (i.e.,
assign the number of occasions using attr(x,"ns")<-NS). This saves
some typing in model specification.
Factors are allowed in ivar and tvar. When a factor is specified,
the contr.treatment coding is used. By default, an intercept
is assumed and the first level of all
factors are dropped from the model (i.e., first levels are the reference levels,
the default R action). However, there are applications where more than one level
will need to be dropped, and the user has control over this via the drop.levels
argument to ivar and tvar. For example,
tvar(x,drop.levels=c(1,2)) drops the first 2 levels of factor x.
tvar(x,drop.levels=length(levels(x))) does the SAS thing and drops
the last level of factor x. If drop.levels is outside the range
[1,length(levels(x))] (e.g., negative or 0), no levels of the factor
are dropped. If no intercept is fitted in the model, this results in the
so-called cell means coding for factors.
Example model specifications: Assume 'age' is a NAN x NS 2-d matrix of ages,
'effort' is a size NS 1-d vector of efforts, and 'sex' is a size NAN 1-d factor
of sex designations ('M' and 'F').
tvar, ivar, print.cjs, residuals.cjs, plot.cjs,
F.cjs.covars, F.cjs.gof, mra.control## Fit CJS model to dipper data, time-varying capture and survivals.
## Method 1 : using factors
data(dipper.histories)
ct <- as.factor( paste("T",1:ncol(dipper.histories), sep=""))
attr(ct,"nan")<-nrow(dipper.histories)
dipper.cjs <- F.cjs.estim( ~tvar(ct,drop=c(1,2)), ~tvar(ct,drop=c(1,6,7)), dipper.histories )
## Method 2 : same thing using 2-d matrices
xy <- F.cjs.covars( nrow(dipper.histories), ncol(dipper.histories) )
for(j in 1:ncol(dipper.histories)){ assign(paste("x",j,sep=""), xy$x[,,j]) } # Extract 2-D matrices of 0s and 1s
dipper.cjs <- F.cjs.estim( ~x3+x4+x5+x6+x7, ~x2+x3+x4+x5, dipper.histories )
## Values in the 1st column of capture covariates do not matter
x3.a <- x3
x3.a[,1] <- 999
dipper.cjs2 <- F.cjs.estim( ~x3.a+x4+x5+x6+x7, ~x2+x3+x4+x5, dipper.histories )
# compare dipper.cjs2 to dipper.cjs
## Values in the last column of survival covariates do not matter
x3.a <- x3
x3.a[,ncol(dipper.histories)] <- 999
dipper.cjs2 <- F.cjs.estim( ~x3+x4+x5+x6+x7, ~x2+x3.a+x4+x5, dipper.histories )
# compare dipper.cjs2 to dipper.cjsRun the code above in your browser using DataLab