########################################################################

# Regression example

########################################################################

model

{

for(i in 1:n)

{

logtime[i] <- log(time[i])

logtime[i] ~ dnorm(mu[i], inv.sigma2)

mu[i] <- b[1] + b[2]*nesting[i] + b[3]*size[i] + b[4]*status[i]

}

b[1] ~ dnorm(0, 0.000001)

b[2] ~ dnorm(0, 0.000001)

b[3] ~ dnorm(0, 0.000001)

b[4] ~ dnorm(0, 0.000001)

inv.sigma2 ~ dgamma(0.0001, 0.0001) # approximates the improper prior

sigma2 <- 1/inv.sigma2

sigma <- sqrt(sigma2)

# prediction of means (mu_p) and new obs (y_p).

for(i in 1:4)

{

mu_p[i] <-b[1]+b[2]*nesting_p[i]+b[3]*size_p[i]+b[4]*status_p[i]

y_p[i] ~ dnorm(mu_p[i], inv.sigma2)

}

# prediction of new obs (ystar) using observed x.

for(i in 1:n)

{

ystar[i] ~ dnorm(mu[i], inv.sigma2)

pright[i] <- step(ystar[i] - logtime[i])

pleft[i] <- step(logtime[i] - ystar[i])

}

}

# Data 1

list(n = 62)

# Data 2

time[] nesting[] size[] status[]

3.030 1.000 0 1

5.464 2.000 0 1

4.098 1.210 0 1

1.681 1.125 0 1

8.850 5.167 0 1

1.493 1.000 0 0

7.692 2.750 0 1

3.846 5.630 0 1

16.667 3.000 0 1

4.219 4.670 0 0

8.130 4.056 0 1

5.000 1.000 0 1

7.299 6.960 0 0

1.000 1.670 0 0

27.027 5.560 0 1

3.106 2.830 0 0

4.000 4.375 0 0

16.129 4.125 0 0

3.484 3.670 0 1

37.037 8.330 0 1

7.299 2.750 0 1

2.525 1.430 0 0

4.132 2.000 0 1

2.000 2.750 0 1

10.000 4.500 0 1

2.667 7.120 0 1

4.587 4.580 0 1

58.824 2.350 0 1

32.258 6.870 1 1

2.571 3.830 1 0

2.160 5.000 1 0

1.000 1.250 1 0

2.967 2.270 1 1

9.524 5.350 1 1

11.111 8.700 1 1

7.299 6.100 1 1

4.000 3.330 1 1

2.381 3.640 1 1

2.611 4.830 1 0

3.257 4.670 1 1

1.701 1.700 1 1

1.795 1.330 1 1

1.198 1.000 1 0

3.185 1.900 1 0

2.273 4.420 1 0

1.111 1.250 1 0

1.000 1.000 1 0

1.000 1.000 1 1

1.230 1.000 1 0

6.061 2.500 1 1

3.175 1.500 1 1

2.000 2.500 1 1

5.076 5.630 1 1

1.934 2.370 1 1

1.493 1.500 1 1

1.000 1.000 1 1

5.102 6.500 1 1

3.003 4.500 1 1

1.898 2.170 1 1

41.667 11.620 1 1

1.000 1.000 1 0

1.000 1.000 1 1

END

# Data 3

nesting_p[] size_p[] status_p[]

4 0 0

4 1 0

4 0 1

4 1 1

END

# Initial values. Initialize other values using "gen inits".

list( b=c(0, 0, 0, 0), inv.sigma2 = 1, y_p=c( 0, 0, 0, 0))