[R-sig-eco] Bayesian analysis of canopy cover/percent data

[Stratford, Jeffrey]

Hi everyone,

I'd like to estimate the means and sd's of canopy cover for 13 sites (9 are
birds and 4 are potential habitats). I'm using a beta distribution to
account for the fact that the data are bound by 0 and 1. For each sample, I
measured the proportion of points that have some sort of vegetation in the
canopy.

I'm trying to adapt the code that works for the other distributions but I'm
getting the error message:

Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) :
  NA/NaN/Inf in 'y'
In addition: There were 33 warnings (use warnings() to see them)
Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) :
  NA/NaN/Inf in 'y'

Below are the R code, the model statement, and the data. I'm using JAGS and
R 3.1.1 in Windows Vista. I suspect the problem is with the model statement.

Any help would be greatly appreciated. I feel like the kid that knows all
the rules of baseball but can neither catch or hit!

Thanks,

Jeff



#########  R  code ##########
frag <- read.csv("f:\\brazil\\TIandFRAG.csv", header=T)
library(R2jags)
library(rjags)
setwd("f://brazil")
site <- frag$site
cover20p <- frag$cover20p/100
N <- length(frag$site)

jags.data <- list("site", "cover20p")
jags.params <- c("median", "test100MF","test100MT","test100fc","test100fa",
"test100gv","test100hm","test100mc", "test100ca","test100ct", "test10MF",
"test10MT", "test10fc","test10fa", "test10gv", "test10hm", "test10mc",
"test10ca", "test10ct",
"test1MF", "test1MT", "test1fc",  "test1fa",  "test1gv", "test1hm",
"test1mc", "test1ca", "test1ct", "t1est1_con","t2est10_con","t3est100_con",
"t4est1_100","t5est1_10","t6est10_100")
#inits1 <- list(a=0, sd=0)
#inits2 <- list(a=100, sd=50)
#jags.inits <- list(inits1, inits2)

jags.inits <- function() {
  list(a=c(0,0,0,0,0,0,0,0,0,0,0,0,0), sd=1)}

jagsfit <- jags(data=jags.data, inits=jags.inits, jags.params,
n.iter=1000000, n.burnin=20000, model.file="fragmodelbeta.txt")

my.coda <- as.mcmc(jagsfit)

summary(my.coda, quantiles=c(0.05, 0.25,0.5,0.75, 0.95))

print(jagsfit, digits=3)

######## MODEL ########################################

model{
for (i in 1:227) {

log(mean[i]) <- a[site[i]]
cover20p[i] ~ dbeta(1,0.5)
}
for (i in 1:13){
a[i] ~ dnorm(0, tau)
median[i] <- exp(a[i])
}
sd ~ dunif(0, 10)
tau <- 1 / (sd*sd) # precision

}

######## DATA #################################################

# normally this would be on a flash drive

structure(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3,
3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
5, 5, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 7, 7,
7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 8, 8, 8, 8, 8, 8,
8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 9, 9, 9, 9, 10, 10, 10, 11, 11,
11, 11, 11, 11, 11, 11, 12, 12, 12, 12, 13, 13, 13, 13, 13, 13,
13, 13, 13, 13, 13, 13, 13, 0.95, 0.8, 0.85, 0.9, 0.35, 1, 1,
1, 0.95, 0.55, 0.9, 0.85, 0.7, 0.65, 0.05, 0.6, 1, 1, 0.85, 0.9,
0, 0.45, 1, 0.7, 0.95, 0.5, 0.95, 0.6, 0.65, 0.7, 0.4, 0.85,
0.6, 0.95, 0.75, 0.9, 0.85, 0.75, 0.7, 0.85, 0.3, 0.7, 0.8, 0.7,
0.75, 0.8, 0.75, 0.95, 0.9, 0.05, 0.85, 0.6, 0.65, 0.5, 0.85,
0.95, 0.85, 0.25, 0.75, 1, 0.65, 0.95, 0.8, 0.9, 0.6, 0.8, 1,
0.2, 0.8, 0.4, 1, 0.95, 0.4, 1, 1, 0.95, 0.45, 0.2, 0.7, 0.95,
0.7, 0.8, 0.5, 0.85, 0.55, 0, 0.25, 0.45, 1, 0.95, 1, 0.9, 0.6,
0.35, 0.95, 0.3, 1, 1, 0.5, 0.4, 0.9, 1, 0.7, 1, 0.9, 1, 0.4,
0.55, 0.8, 0.7, 1, 0, 0.8, 0, 0.7, 0.5, 0.8, 0.75, 0, 0.45, 0.1,
0, 0.4, 0.55, 0.4, 1, 0.9, 0.9, 0.15, 0.55, 0.35, 0.9, 0.65,
0.25, 1, 0.85, 1, 0.95, 0.7, 0.5, 0.7, 0.2, 0.95, 1, 1, 0.25,
0.85, 0.5, 0.8, 0.75, 0.85, 0.7, 0.95, 0.05, 0.65, 0.65, 1, 1,
1, 0.65, 0.4, 0.6, 0.9, 0.85, 0.75, 0.5, 0.65, 1, 0.65, 0.55,
0.75, 0.4, 0.9, 0.35, 1, 1, 0.4, 0.5, 0.8, 0.95, 0.95, 0.55,
0.7, 0.85, 0.8, 0.8, 0.65, 1, 0.6, 0.5, 1, 0.8, 1, 0.45, 1, 1,
0.8, 0.85, 1, 1, 1, 1, 0.5, 0.6, 0.15, 0.75, 0.6, 0.1, 0.05,
0, 1, 0.6, 0.1, 0.35, 0.9, 0.9, 0.95, 0.95, 0.9, 0.55, 0.65,
0.9, 0.4, 1, 0.65, 0.5, 0.8), .Dim = c(227L, 2L), .Dimnames = list(
    NULL, c("site2", "cover20p2")))



-- 
********************************************************
Jeffrey A. Stratford, PhD
Department of Biology and Health Sciences
Wilkes University, PA 18766 USA
570-332-2942
http://web.wilkes.edu/jeffrey.stratford/
********************************************************

	[[alternative HTML version deleted]]