Estimates of error - Package marked

[kbearden]

Good evening,

I am working on a multi-state model with 8 states and 1 environmental covariate. I am working in package marked and the main issue that I am running into is that I am not getting any estimates of error for the parameters (please see below the code for a model that I have run). I would greatly appreciate any help and or direction. Please let me know if more code or context is needed.

Code: Select all

#### mod.ss.pt.psiorXs ####
mod.ss.pt.psiorXs <- crm(pdfh1.proc, ddl = pdfh1.ddl, model = "Mscjs",
                        model.parameters = list(p=list(formula = ~ time),
                                                S=list(formula = ~ stratum),
                                                Psi=list(formula = ~ or*stratum)),
                        burnin = 100, iter = 1000, use.admb=TRUE)
mod.ss.pt.psiorXs.real <- predict(mod.ss.pt.psiorXs, ddl = pdfh1.ddl)
mod.ss.pt.psiorXs.real

I should note that "or" stands for open river and it has been added to my design data.

Thank you,
Katie

[jlaake]

set hessian=TRUE in the crm function call. Unlike MARK, the default doesn't compute the hessian by default to save computation time. MARK uses the hessian to count estimable parameters which is useful but not infallible. I chose to not estimate the hessian to speed up model selection and assume the user is able to ascertain the correct number of parameters - for better or worse. Once you have selected a model, then you can set hessian=TRUE to get error estimates. Note that burnin is only used for the one Bayesian model probitcjs.

[kbearden]

Thank you jlaake for the response / help. I have updated my code (below) but am still not getting any estimates of error.

Code: Select all

#### mod.ss.pt.psiorXs2 - with Hessian = TRUE and NO burnin ####
mod.ss.pt.psiorXs2 <- crm(pdfh1.proc, ddl = pdfh1.ddl, model = "Mscjs",
                          model.parameters = list(p=list(formula = ~ time),
                                                  S=list(formula = ~ stratum),
                                                  Psi=list(formula = ~ or*stratum)),
                          iter = 1000, use.admb=TRUE,
                          hessian = TRUE)
mod.ss.pt.psiorXs2.real <- predict(mod.ss.pt.psiorXs2, ddl = pdfh1.ddl)
mod.ss.pt.psiorXs2.real

Here is a small portion of the output I am getting when I run the above code.

Code: Select all

$Psi
    stratum tostratum or occ     estimate
1         A         A  0  45 1.000000e+00
2         A         B  1  45 1.617269e-03
3         A         C  2  45 2.700116e-03
4         A         D  0  45 9.686836e-04
5         A         E  1  45 1.617269e-03
6         A         F  0  45 9.686836e-04
7         A         G  0  45 9.686836e-04
8         A         H  0  45 9.686836e-04
9         B         A 12  45 1.477142e-02
10        B         B 12  45 1.000000e+00

Any additional insight is greatly appreciated!

Thank you,
Katie

[jlaake]

Katie

You had said " I am not getting any estimates of error for the parameters". I assumed that you meant the beta parameters of the model. Adding hessian=TRUE should have produced that. If what you want is error estimates from predict you need to add se=TRUE in that call. This is described in the help ?predict.crm.
Also there is an example("mstrata") that shows various ways that mscjs can be fitted. I no longer have admb on my computer because TMB is a much better approach. See use.tmb=TRUE or you can set both to FALSE and it uses std optimization routine. Please look through the help files to see what arguments you can set to control behavior of the functions.

--jeff

[kbearden]

Hi Jeff,

I was able to get what I was looking for by adding se = TRUE in the predict function. Thank you very much for your help and guidance!

Best,
Katie