profile confidence intervals in RMark

posts related to the RMark library, which may not be of general interest to users of 'classic' MARK

profile confidence intervals in RMark

Postby TimG » Mon Jan 08, 2018 12:24 pm

I am using the Huggins CRDMS model in RMark. I have a cumbersome dataset - 22 primary periods, 3 states (1 of which is unobserable with p=0), sex as group, 'age model' based on age at marking - so working in the MARK GUI is not really practical.

I can fit a simple model [S~1, p=c~1, Psi~stratum], but get the following error when I set profile.int=TRUE for the same model:
Code: Select all
Error in wtable[i, cols] = real[pims[i, cols]] :
 number of items to replace is not a multiple of replacement length

I get this message regardless of whether I specify mlogit0=TRUE or options="SIMANNEAL". Some parameters are fixed (by adding the ddl$Psi$fix field to the design data), so not sure if that is affecting things.

Can anyone provide insight on this?
TimG
 
Posts: 8
Joined: Mon Aug 31, 2015 12:37 pm

Re: profile confidence intervals in RMark

Postby jlaake » Mon Jan 08, 2018 1:21 pm

Look at the .out file from Mark to see what it contains.
jlaake
 
Posts: 1020
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA

Re: profile confidence intervals in RMark

Postby TimG » Mon Jan 08, 2018 4:55 pm

The last line in the .out file is "Error -- length of command exceeded MAXCHR/2 characters".
This occurs after the PIMs for S, Psi, p, and c are specified (and after info for initial values if I include that option in the mark() call).

It looks like Mark created the INP file ok, but nothing is being fit since I get that message in R almost immediately and no RES or VCV files are created.
TimG
 
Posts: 8
Joined: Mon Aug 31, 2015 12:37 pm

Re: profile confidence intervals in RMark

Postby jlaake » Tue Jan 09, 2018 3:08 pm

Tim sent me his files offlist. There were 16,000+ real parameters in the model and what is passed to MARK is the indices for each if you set profile.int=TRUE. This exceeded the maximum line length that MARK would accept. Instead, specify a vector of indices for the real parameters that you want profile confidence intervals. Just use a vector with far fewer than 16,000 indices.

--jeff
jlaake
 
Posts: 1020
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA

Re: profile confidence intervals in RMark

Postby TimG » Fri Jan 12, 2018 11:42 am

Thanks, Jeff - that did it!

Extracting the parameter indices was a bit tricky for more complex models, so here’s what I came up with, in case it might be helpful for others:
-Run the model with “profile.int=FALSE”. Here, I’ve named this model “m1”.
-“unique(m1$results$real)” will return the unique, real parameter estimates, but the “get.real()” function returns the par.index field you’ll need, so join these tables:
Code: Select all
merge(x=unique(m1$results$real), y=get.real(m1,"Psi",se=TRUE), by="row.names")

Note that the “unique()” function may not return some parameters if they were estimated to be the same value as another parameter (eg. both were at a boundary, or had the same covariate value in a covariate model), so you may have to double-check that all parameters of interest are included.
-Repeat this join for each parameter type, as needed (in my case: S, p, and Psi).
-Re-run the model, specifying the vector of par.index values for profile CIs: eg. “profile.int=c(1,2,3)”.
Be patient as the model re-runs!
TimG
 
Posts: 8
Joined: Mon Aug 31, 2015 12:37 pm


Return to RMark

Who is online

Users browsing this forum: No registered users and 1 guest