inestimable parameters differ between Mark and RMark

[aber_r16]

HI All, before I begin I have read the forum post by Jeff regarding differences in results between Mark and RMark...same data (checked), think the design matrix are the same (this could be the problem), and checked the link function (logit in both cases). I can't confirm that the design matrix are the same because I built the MARK model with PIMS (should be the same, but perhaps I will learn Im wrong:)

Im doing a TSM CJS model. The estimable parameters are the same when comparing models - its just the inestimable parameters that differ. As I know these estimates are not usable I am tempted to forge ahead without worrying about this, but its bugging me that the results aren't exactly the same. Is there some theoretical explanation or am I wrong about meeting the three conditions above (same data, same design matrix, and same links).

Cheers! Robin

[jlaake]

Robin-

You didn't provide enough detail for us to evaluate but if you created the model with PIMS in MARK interface then it would have used the sin link by default whereas RMark would use logit link. Did you force the logit link in MARK? You'll have to provide a whole lot more details before anyone can be of any help. Did you try and export to MARK with export.MARK and run the models from the MARK interface and compare the models side by side? I don't know what you mean by inestimable parameters but if it is something like Phi(t)p(t) product at last occasion then the individual values can differ but the product should be the same. Were the likelihoods the same? More details.

--jeff

[aber_r16]

Hi Jeff, thanks for the INSTANT reply! Its so helpeful to get a timely response! Sorry for the lack of detail but you did manage to answer my question.

To backup - yes, I forced MARK to use the logit link so I knew I was comparing apples to apples.

By inestimable parameters I did mean the Phi and p parameters for the last occasion. However, when I looked at each model (MARK) vs (RMark) and multiplied Phi*p I did indeed get the same results so am satisfied I have indeed run the same model via each method. Is there an explanation of why the "confounded" (is that what you would call them) parameters would be different in Mark vs RMark?

Thanks! Robin

[jlaake]

Technically any pair of values that give the same product will have the same likelihood value. So you could end up with different values by simply changing the starting values for the beta parameters. Since mark.exe is used in both the MARK and RMark interfaces, there must be a difference in the design matrix (DM) that is being used. Given that you used PIMS in the MARK interface to build the model, the DM will be the identity matrix. I don't know what you did with RMark but unless you used a formula that removed the intercept the DM will not be an identity matrix and that could explain the difference. If the model we are talking about is Phi(time)p(time), use the formula ~-1+time for both Phi and p and then RMark will create an identity DM and they should be the same. Also if you do that then you could use link="sin" and then compare that to what MARK gives with the PIM and sin link. In theory, any DM that specifies the same real parameter values should yield the same likelihood and beta estimates. But this is a numerical procedure and small differences can occur. In this case, there is no unique solution for those parameters so any pair that produces the same product of the real values is a valid solution.

--jeff

[aber_r16]

Hi Jeff,

I tinkered with this problem a little bit more and determined that you are absolutely correct - of course! The design matrix were different, and its all do do with how MARK vs RMark handles the intercept. MARK codes it last, RMark codes it first.

Is there a simple way to change this in RMark?

Thanks for your patience with these beginner questions!!

-Robin

[jlaake]

Correct. RMark is using the R convention (of course) setting the first level to the intercept and MARK is using the SAS convention of using the last level. Use the relevel function in R to set the intercept to the last level.

Something like

ddl$p$time=relevel(ddl$p$time,"k")

where k is the value of the last time. You can set any factor level as the intercept. I think it is best to use the level with the most data. Using the level with very sparse or no data as the intercept will get you a bunch of nonsense std errors because the intercept is used for each parameter.

--jeff