rbchan / unmarked

R package for hierarchical models in ecological research
https://rbchan.github.io/unmarked/
37 stars 25 forks source link

error in ranef() for multmixOpen objects #204

Open dslramsey opened 3 years ago

dslramsey commented 3 years ago

Missing values in y for multmixOpen() models causes ranef() to fail. A reproducible example follows.

set.seed(123)
lambda=4; gamma=0.5; omega=0.8; p=0.5
M <- 100; T <- 5
y <- array(NA, c(M, 3, T))
N <- matrix(NA, M, T)
S <- G <- matrix(NA, M, T-1)

for(i in 1:M) {
  N[i,1] <- rpois(1, lambda)
  y[i,1,1] <- rbinom(1, N[i,1], p)    # Observe some
  Nleft1 <- N[i,1] - y[i,1,1]         # Remove them
  y[i,2,1] <- rbinom(1, Nleft1, p)   # ...
  Nleft2 <- Nleft1 - y[i,2,1]
  y[i,3,1] <- rbinom(1, Nleft2, p)

  for(t in 1:(T-1)) {
    S[i,t] <- rbinom(1, N[i,t], omega)
    G[i,t] <- rpois(1, gamma)
    N[i,t+1] <- S[i,t] + G[i,t]
    y[i,1,t+1] <- rbinom(1, N[i,t+1], p)    # Observe some
    Nleft1 <- N[i,t+1] - y[i,1,t+1]         # Remove them
    y[i,2,t+1] <- rbinom(1, Nleft1, p)   # ...
    Nleft2 <- Nleft1 - y[i,2,t+1]
    y[i,3,t+1] <- rbinom(1, Nleft2, p)
  }
}
y=matrix(y, M)

y[1,1]<- NA   # Insert NA (could be anywhere)

#Create some random covariate data
sc <- data.frame(x1=rnorm(100))

## Not run: 
#Create unmarked frame
umf <- unmarkedFrameMMO(y=y, numPrimary=5, siteCovs=sc, type="removal")

#Fit model
fit <- multmixOpen(~x1, ~1, ~1, ~1, K=30, data=umf)

# Empirical Bayes estimates of abundance for each site / year
re <- ranef(fit)

Error in if (any(x < 0)) stop("'x' must be non-negative") : 
missing value where TRUE/FALSE needed 

The issue appears to be in the call to dmultinom() where NA's are trapped for capture probability but not for y.

kenkellner commented 3 years ago

This specific situation should be fixed in 7fc9ca5e534ccb1d423a56edaf9d2ea4f8626b45. I'll leave the issue open because I haven't addressed the situation where an entire year of observations are NA, which gets complicated for this model.

dslramsey commented 3 years ago

Sorry for the radio silence. Just got around to testing the latest changes and the issue has been now resolved for me.