Occupancy GOF discrepancy between Bayesian and Frequentist

[Luca Petroni]

Dear all,

I am running a single-season single-species Bayesian occupancy model on systematic camera trapping data calling JAGS from jagsUI in R.
I planned of testing GOF by using the code for the "closed part of the model" from chapter 4.8.2 in AHM2 (Kery and Royle, 2021). Basically, it is a calculation of both a Chi-squared and a Freeman-Tukey discrepancy stats. 
The problem is, I obtain a very bad fit and overdispersion for the NULL MODEL, despite having a quite robust dataset (33 sites simultaneously active for 3 months, naive occupancy approx. 0.6). Furthermore, and this is the strangest part, if I try fitting the same null model in unmarked, both the MacKenzie-Bailey test and the same Chi-squared and FT statistics (code from Andrade-Ponce et al. 2022) result in a good fit and no overdispersion at all.

Is there a plausible reason why I obtain such different results in GOF testing between Bayesian and Frequentist methods? Am I missing something?

Thank you in advance.

Luca

[Matthijs Hollanders]

Can you share the code you’ve used in JAGS?

--
*** Three hierarchical modeling email lists ***
(1) unmarked: for questions specific to the R package unmarked
(2) SCR: for design and Bayesian or non-bayesian analysis of spatial capture-recapture
(3) HMecology (this list): for everything else, especially material covered in the books by Royle & Dorazio (2008), Kéry & Schaub (2012), Kéry & Royle (2016, 2021) and Schaub & Kéry (2022)
---
You received this message because you are subscribed to the Google Groups "hmecology: Hierarchical Modeling in Ecology" group.
To unsubscribe from this group and stop receiving emails from it, send an email to hmecology+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/hmecology/29824749-62d7-480a-b5c5-d6eedf000caen%40googlegroups.com.

--

Dr. Matthijs Hollanders

Statistical Consultant – Quantecol

Postdoctoral Research Fellow – College of Science | Australian National University

[Luca Petroni]

Hi Dr. Hollanders,

here is the code, almost the same as in the book.

# Draw a replicate data set under the fitted model
for (i in 1:nsites){
    for (j in 1:nsurveys){
      yrep[i,j] ~ dbern(z[i] * p)
  }
}

# Compute detection frequencies for observed and replicated data
for (i in 1:nsites){
    # Det. frequencies for observed and replicated data
    detfreq[i] <- sum(y[i,])
    detfreqrep[i] <- sum(yrep[i,])
    # Expected detection frequencies under the model
    for (j in 1:nsurveys){
      tmp[i,j] <- z[i] * p
    }
    E[i] <- sum(tmp[i,])     # Expected number of detections
    # Chi-square and Freeman-Tukey discrepancy measures
    # ..... for actual data set
    x2[i] <- pow((detfreq[i] - E[i]),2) / (E[i]+0.00000001)
    ft[i] <- pow((sqrt(detfreq[i]) - sqrt(E[i])),2)
    # ..... for replicated data set
    x2rep[i] <- pow((detfreqrep[i] - E[i]),2) / (E[i]+0.00000001)
    ftrep[i] <- pow((sqrt(detfreqrep[i]) - sqrt(E[i])),2)
  }

# Add up Chi-square and FT discrepancies and compute fit stat ratio (closed part)
Chi2 <- sum(x2)
FT <- sum(ft)
Chi2rep <- sum(x2rep)
FTrep <- sum(ftrep)
Chi2.chat <- Chi2 / Chi2rep
FT.chat <- FT / FTrep
Chi2.bpv <- step(Chi2rep-Chi2)
FT.bpv <- step(Chi2rep-Chi2)

Thanks for your response,

Luca