<html>
  <head>
    <meta content="text/html; charset=ISO-8859-1"
      http-equiv="Content-Type">
  </head>
  <body text="#000000" bgcolor="#FFFFFF">
    <div class="moz-cite-prefix">On 12-11-07 11:32 AM, Ben Bolker wrote:<br>
      <br>
      Well I don't recall getting personal. I was referring to<br>
      what I feel are problems with the design of R.<br>
      <br>
      Anyway it apears that you don't seem to have the right model here.<br>
      <br>
      This should be pretty simple.  One has 192 random effects.  all I
      did<br>
      was to use 3 different parameters for the standard deviations of
      the random<br>
      effects.  sd1 for obs 1-64  sd2 for obs 65-128 and sd3 for obs
      129-192.<br>
      You seem to have something completely different.  the par file is<br>
      <br>
      # Number of parameters = 8  Objective function value = 285.193 
      Maximum gradient component = 6.97992e-05<br>
      # pz:<br>
      0.000100000000000<br>
      # beta:<br>
       -1.68901008026 -0.817412602915 -4.06986974939 6.84143453412<br>
      # tmpL:<br>
       -0.114257455234 -0.256577632335 -0.00360094954491<br>
      # tmpL1:<br>
       0.000100000000000 0.000100000000000 0.000100000000000<br>
      # log_alpha:<br>
      1.36928124319<br>
      <br>
      so the fit is much worse.<br>
      <br>
      <br>
      <br>
      <br>
      <br>
    </div>
    <blockquote cite="mid:509AB762.9050801@gmail.com" type="cite">
      <pre wrap="">On 12-11-07 01:17 PM, dave fournier wrote:
</pre>
      <blockquote type="cite" style="color: #000000;">
        <pre wrap=""><span class="moz-txt-citetags">> </span>But that is not right.  The model is to assume that the
<span class="moz-txt-citetags">> </span>random effects in different time periods are independent and that their
<span class="moz-txt-citetags">> </span>variances are different.  You can do this by figuring out what the
<span class="moz-txt-citetags">> </span>correct covariance is for the current design matrix or by
<span class="moz-txt-citetags">> </span>writing the model properly which is trivial to do in ADMB.  As usual it
<span class="moz-txt-citetags">> </span>is the R interface and these silly design matrices which introduces the
<span class="moz-txt-citetags">> </span>problem.
<span class="moz-txt-citetags">> </span>_______________________________________________
<span class="moz-txt-citetags">> </span>Users mailing list
<span class="moz-txt-citetags">> </span><a moz-do-not-send="true" class="moz-txt-link-abbreviated" href="mailto:Users@admb-project.org">Users@admb-project.org</a>
<span class="moz-txt-citetags">> </span><a moz-do-not-send="true" class="moz-txt-link-freetext" href="http://lists.admb-project.org/mailman/listinfo/users">http://lists.admb-project.org/mailman/listinfo/users</a>
</pre>
      </blockquote>
      <pre wrap="">
  Just to continue this a little bit:

If I do this:

l2 <- read.csv("bigelow_dat.txt")

## unnecessary picture
library(ggplot2)
ggplot(l2,aes(factor(time),snag))+
    geom_dotplot(binaxis="y",stackdir="center")+
    facet_wrap(~Ecozone+time,scale="free_x")

with(l2,table(table(plot))) ## check exp. design: 3 samples per plot
with(l2,table(table(plot,time))) ## check exp. design: 1 sample per
plot/time
l2$ftime <- factor(l2$time) ## factor-ize time

library(glmmADMB)
## note zero-inflation is FALSE/off throughout (default)
g1 <- glmmadmb(formula = snag ~ ftime + Ecozone + (1 | plot), data = l2,
family = "nbinom1")
VarCorr(g1)
log(sqrt(VarCorr(g1)[[1]]))  ## log std. dev. =0.4315505  ## matches DF
logLik(g1)
##  'log Lik.' -255.758 (df=6)  ## matches DF
coef(g1)
## (Intercept)      ftime2      ftime3    EcozoneW
##   -0.7635323  -0.5884655  -0.6037133   1.2493021

## ?? The fixed effects don't match DF results: did he orthogonalize the
design matrix ... ?
## never mind

## besides being the wrong model, this fails on my system (32-bit Linux)
try(g2 <- glmmadmb(formula = snag ~ ftime + Ecozone + (ftime | plot),
data = l2, family = "nbinom1"))

## this works and is (I think) the right model ...
g3 <- glmmadmb(formula = snag ~ ftime + Ecozone + (ftime-1 | plot), data
= l2, family = "nbinom1",
               save.dir="glmmadmb_tmp")
## the Z matrix looks like I think it should (an indicator variable for
time period)
VarCorr(g3)

## Random effects:
## Structure: Diagonal matrix
## Group=plot
##       Variance StdDev
## time1   0.7957 0.8920
## time2   0.5986 0.7737
## time3   0.9928 0.9964

log(sqrt(diag(VarCorr(g3)[[1]])))
##      ftime1      ftime2      ftime3
## -0.11426024 -0.25658084 -0.00360295

## This doesn't match DF: he gets the log-sd of the variances as

## tmpL:
##  0.352458239856 0.489012745656 0.496221176099

## i.e. increasing from period 1 to periods 2 and 3

For the overdispersion model (which I haven't run yet -- glmmADMB is not
flexible enough) -- DF gets

# tmpL:
 0.436593520224
# tmpL1:
 0.000100000000000
# log_alpha:
 0.826338727816 0.0677913813376 0.657156858701

i.e. the log-sd among plots (tmpL) is about the same as before, and
log_alpha is the log of the NB scale parameter.  This suggests <b class="moz-txt-star"><span class="moz-txt-tag">*</span>lower<span class="moz-txt-tag">*</span></b>
variance (at least when scaled by the mean) in the second period,
<b class="moz-txt-star"><span class="moz-txt-tag">*</span>higher<span class="moz-txt-tag">*</span></b> variance in the first period -- which seems (??) to match the
(ftime-1|plot) glmmADMB results rather than the DF results (??)

  I'd rather be arguing over what the data say than over whose tools are
better
or worse, or who's an idiot.  I'm happy to concede that I'm an idiot, so we
can get that over right at the beginning ...

  cheers
</pre>
    </blockquote>
    <br>
  </body>
</html>