[R-sig-ME] glm module in

[Ben Bolker]

On 09/26/2011 06:13 AM, Martyn Plummer wrote:
> On Tue, 2011-09-20 at 08:50 -0400, Ben Bolker wrote:
>> On 09/20/2011 05:34 AM, Martyn Plummer wrote:
>>> Hi Ben,
>>>
>>> There are two reasons that your model is not recognized by the the glm
>>> module. I'll repeat the essential parts of the model:
>>>
>>> SibNeg[i] ~ dpois(mu[i])
>>> mu[i] <- lambda[i] * z[i] + 0.00001
>>> log(lambda[i]) <- offset[i] + alpha + inprod(X[i,],beta) + a[nest[i]]
>>>
>>> Firstly, this is not a log-link GLM due to the perturbation you
>>> have added to mu[i].  Secondly, and this is the real problem, even
>>> without the perturbation JAGS still would not recognize it. 
>>>
>>> Normally you can work around this problem by adding log(z) as an offset,
>>> but in your case this will not work since z may be zero.
>>>
>>> I am aware of this issue: it is a commonly occurring set-up in
>>> epidemiological applications, where lambda would a be in incidence rate
>>> and z would be the person-time at risk. I'll put this on the TODO list.
>>>
>>> Martyn
>>
>>   Thanks, Martyn, for the quick response.  Nice to know I wasn't missing
>> something too obvious.  It's not critical for this application (I can
>> afford to let the model run for 10 minutes), but I wonder if a
>> workaround would be to use something like
>>
>>    logzzvals <- c(-10,0)
>>    lzcat[i] ~ dbern(prob)
>>    log(lambda[i]) <- offset[i] + logzvals[lzcat[i]+1] + ...
>>
>> ?
> 
> I'm afraid not.   An offset normally works in this situation, but if
> z[i] == 0 then SibNeg[i] is completely uninformative about the
> parameters alpha, beta, and a. 

  Yes, but here I designed it so z[i] was never zero, but was either
exp(-10) or 1.  (Again, I suppose I could just try this and see if it
worked.)


  thanks,
     Ben