RMark Huggins with c=p+offset

[jlaufenb]

I assume that you are attempting to build a model with an additive behavioral effect by listing paramter-specific submodels for p and c. If so, then try the following insted:

pPlusc<-list(formula=~1+c, share=TRUE)

#with only pPlusc in your parameter list to code for p and c; i.e., model.paramters=list(p=pPlusc)

I use this with the FullHet datatype and believe is should work with the HugFullHet datatype, too.

[timkdavies]

I changed the code to what you suggested, unfortunately the same ridiculous estimates, se and ci still came out of the model. It seems that every time c is included in the formula the estimates go a bit hay-wire.

The top 2 models both have c as a variable -
model npar AICc DeltaAICc weight
7 pi(~1)p(~mixture + time + c)c() 10 561.1849 0.000000 9.769305e-01
11 pi(~1)p(~time + c)c() 9 568.8016 7.616738 2.167255e-02
14 pi(~1)p(~time + mixture)c() 9 574.2852 13.100298 1.396916e-03

I'm quite confused as to why c would have such a huge impact on the estimates?

These are the top 2 rows of estimates from the top 4 models (lowest AIC values)

Model 7-
Group N-hat Standard Error Lower Upper
--------- -------------- -------------- -------------- --------------
1 0.1093380E+012 0.3018018E+012 0.6158692E+010 0.1941128E+013
2 0.5466902E+011 0.1557737E+012 0.2966273E+010 0.1007561E+013

Model 11-
Group N-hat Standard Error Lower Upper
--------- -------------- -------------- -------------- --------------
1 7.5080574 16.650958 2.2789217 110.77138
2 3.7540287 8.6303505 1.1337779 57.696001

Model 14-
Group N-hat Standard Error Lower Upper
--------- -------------- -------------- -------------- --------------
1 3.6017476 3.9123595 2.1043912 26.576737
2 1.8019940 2.1350468 1.0471518 14.640914


Model 15- pi(~1)p(~time)c()
Group N-hat Standard Error Lower Upper
--------- -------------- -------------- -------------- --------------
1 2.0750044 0.2798602 2.0029899 3.8815591
2 1.0375496 0.1977010 1.0010348 2.3625755

Should I just be rejecting these models as they are not realistic? Or is there something going on with my code or data?

Cheers!

[timkdavies]

I changed the code to what you suggested, unfortunately the same ridiculous estimates, se and ci still came out of the model. It seems that every time c is included in the formula the estimates go a bit hay-wire.

The top 2 models both have c as a variable -

Code: Select all

                             model npar     AICc  DeltaAICc       weight
7  pi(~1)p(~mixture + time + c)c()   10 561.1849   0.000000 9.769305e-01
11           pi(~1)p(~time + c)c()    9 568.8016   7.616738 2.167255e-02
14     pi(~1)p(~time + mixture)c()    9 574.2852  13.100298 1.396916e-03
15               pi(~1)p(~time)c()    8 597.6915  36.506642 1.154907e-08

I'm quite confused as to why c would have such a huge impact on the estimates?

These are the top 2 rows of estimates from the top 4 models (lowest AIC values)

Code: Select all

Model 7- Group     N-hat        Standard Error      Lower           Upper
 ---------  --------------  --------------  --------------  --------------
       1    0.1093380E+012  0.3018018E+012  0.6158692E+010  0.1941128E+013
       2    0.5466902E+011  0.1557737E+012  0.2966273E+010  0.1007561E+013

Model 11- Group     N-hat        Standard Error      Lower           Upper
 ---------  --------------  --------------  --------------  --------------
       1    7.5080574       16.650958       2.2789217       110.77138     
       2    3.7540287       8.6303505       1.1337779       57.696001       
 
Model 14- Group     N-hat        Standard Error      Lower           Upper
 ---------  --------------  --------------  --------------  --------------
       1    3.6017476       3.9123595       2.1043912       26.576737   
       2    1.8019940       2.1350468       1.0471518       14.640914
 
Model 15- Group     N-hat        Standard Error      Lower           Upper
 ---------  --------------  --------------  --------------  --------------
       1    2.0750044       0.2798602       2.0029899       3.8815591     
       2    1.0375496       0.1977010       1.0010348       2.3625755 Model 7-
 Group     N-hat                     Standard Error      Lower               Upper
 ---------  --------------  --------------  --------------  --------------
       1    0.1093380E+012    0.3018018E+012     0.6158692E+010  0.1941128E+013
       2    0.5466902E+011    0.1557737E+012     0.2966273E+010  0.1007561E+013

Should I just be rejecting these models as they are not realistic? Or is there something going on with my code or data?

Cheers!

[vankatwy]

Hi Tim,
I am fairly new to RMark notation so take my suggestion with a grain of salt. I noticed that c is specified in the p brackets (i.e., pi(~1)p(~mixture + time + c)c()), but if you want them both to be constant and not equal to each other are you not looking for pi(~1)p(~1)c(~1)?
From what I understand the pi(~1)p(~1)c(~1) model is the default for the all closed model types so you could simply use the default model:

model = mark(data, model="FullHet", model.parameters=list())

This approach works for me when creating the pi(~1)p(~1)c(~1) model in Huggins data type.

--Kristin