[ADMB Users] Integrating R and ADMB to generalize TPL files
dave fournier
davef at otter-rsch.com
Wed Nov 7 09:49:27 PST 2012
On 12-11-07 09:29 AM, Ben Bolker wrote:
> On 12-11-07 11:58 AM, dave fournier wrote:
>> Here is an example of why I dislike the R approach
>> so much.
>>
>> https://stat.ethz.ch/pipermail/r-sig-mixed-models/2012q4/019349.html
>>
>> I was curious to find out why introducing to additional parameters for the
>> variance of the random effects would produce no improvement in the
>> model. So first I wanted to run glmmadmb as seth had done.
>>
>> I put his data in a file and read it into R.
>>
>> Then I ran his model by grabbing his call to glmmadmb.
>> I was surprised to see that I got a different result from Seth.
>>
>> Examination of the design matrix and I realized that one of the columns of
>> the data needed to be changed to a factor. If I had not got the
>> different result I would have thought that the analysis was correct.
>> That was the first gotcha. With R you never really know what
>> you are doing. This was for the simpler model with paramerter for
>> the variance of the RE's.
> I did see that, and I agree there is a price to pay for convenience.
It is not convenient to effortless get to the wrong answer,
just a bit of a waste of time.
>> Then I used Seth's R script to run the more complicated model.
>> He was right that the likelihood did not change. Strange!
>>
>> So I examined the design matrix for the RE's.
>> The problem was obvious. rather than a different
>> parameter for each time period the design matrix had
>> an overall mean and differences from it. While this
>> would be the same for fixed effects it is not equivalent
>> for random effects because it assumes a different covariance
>> structure.
>>
>> I hacked together the correct model by hand and ran it on his data.
> I agree that the model Seth was specifying was that there is variation
> among trees in their time-1 response; variation in the difference
> between time 1 and time 2; and variation in the difference between time
> 1 and time 3 (unless I'm wrong, which has been known to happen).
All I know is that he said
The data suggest that variance may differ
sharply before and after treatment (factor variable 'time', values 1/2/3),
I don't see anything about differences just that the variance should
depend on time.
>
> It sounds like you fitted a model with variance around the time-1
> response, variance around the time-2 response, and variance around the
> time-3 response. I *think* one could specify that via
> random=~time-1|subject, but I haven't tried it.
>
> I agree that unless one specifies corStruct="full" these are not
> equivalent models.
>
> Do we know which of these two models answers the question Seth wanted
> to ask?
>
>
>> Results are:
>>
>> # Number of parameters = 8 Objective function value = 255.508 Maximum
>> gradient component = 1.44753e-06
>> # pz:
>> 0.000100000000000
>> # beta:
>> -10.7785676169 -2.55267985235 -4.72926437596 6.91095066201
>> # tmpL:
>> 0.352458239856 0.489012745656 0.496221176099
>> # tmpL1:
>> 0.000100000000000
>> # log_alpha:
>> 0.591420283490
>> the original fit is
>>
>> # Number of parameters = 6 Objective function value = 255.758
>> Maximum gradient component = 2.23893e-05
>> # pz:
>> 0.000100000000000
>> # beta:
>> -10.4061157821 -1.87212686310 -3.41511253652 8.12797421786
>> # tmpL:
>> 0.431551628285
>> # tmpL1:
>> 0.000100000000000
>> # log_alpha:
>> 0.593869268695
>>
>> So it appears that the variance is larger in the second and third time
>> period, but not
>> significantly so.
>>
>> I welcome your snide comments.
>>
>>
>>
>> _______________________________________________
>> Users mailing list
>> Users at admb-project.org
>> http://lists.admb-project.org/mailman/listinfo/users
More information about the Users
mailing list