[R] bugs() ignores my inits

toby909 at gmail.com toby909 at gmail.com
Sat Oct 27 00:01:23 CEST 2007


Hi All

I can specify whatever inits, it has no effect on the estimation. I am 
replicating a textbook example. The result is completely trash, having estimates 
of -58.7 (sd=59.3), where it should be closer to an ml estimate of 0.585 (SE=0.063).

The two chains within one run are different, but with different inits for 
different runs, I get exactly the same chains, and I mean exactly.

If I set strong (high-precision) priors using the ML estimates, I get the result 
I want. I dont want to set strong priors, but I want to set reasonable inits; 
ones that are not taken by sampling from the prior, but ones that I! specify.

Thanks Toby



data <- list("n", "n3", "y", "person", "occ")

model <- function() {
...
}
write.model(model, "example.txt")

inits <- list(list(m1=0.585, m2=0.718, m3=0.672, m4=0.639),list(m1=0.585, 
m2=0.718, m3=0.672, m4=0.639))
parameters <- c("m1", "m2", "m3", "m4", "sig21", "sig22")

bugs1 <- bugs(data, inits, parameters, "example.txt", n.chains=2, n.iter=15000, 
clearWD=TRUE, debug=1, DIC=0)






display(log)
check(J:/project/ps/code/example.txt)
model is syntactically correct
data(J:/project/ps/code/data.txt)
data loaded
compile(2)
model compiled
inits(1,J:/project/ps/code/inits1.txt)
this chain contains uninitialized variables
inits(2,J:/project/ps/code/inits2.txt)
this chain contains uninitialized variables
gen.inits()
initial values generated, model initialized
thin.updater(15)
update(500)
set(m1)
set(m2)
set(m3)
set(m4)
set(sig21)
set(sig22)
update(500)
coda(*,J:/project/ps/code/coda)
stats(*)

Node statistics
	 node	 mean	 sd	 MC error	2.5%	median	97.5%	start	sample
	m1	-58.66	59.27	9.001	-120.8	0.4013	0.7102	501	1000
	m2	-58.53	59.27	9.001	-120.7	0.5417	0.8314	501	1000
	m3	-58.58	59.27	9.0	-120.7	0.4783	0.7871	501	1000
	m4	-58.6	59.27	9.0	-120.7	0.3555	0.7695	501	1000
	sig21	0.07409	0.01034	3.039E-4	0.05636	0.07319	0.09558	501	1000
	sig22	7206.0	7500.0	1097.0	0.07909	8436.0	19810.0	501	1000
history(*,J:/project/ps/code/history.odc)

History

save(J:/project/ps/code/log.odc)
save(J:/project/ps/code/log.txt)



More information about the R-help mailing list