Dear list, I have troubles with using Bayesian logistic model with count data in bayesglm. If I consider the following artificial data set, with a response y and a covariate x, in a form of row data (a) and count (b) :
a<-data.frame(y=c(rep(1,10),rep(0,6)),x=c(rep(5,6),rep(4,4),rep(5,1),rep(4,5))) a$un<-1 b<-aggregate(a$un,list(y=a$y,x=a$x),sum) names(b)[3]<-"w" A binomial model with glm running on a and b give the same results : (M1=glm(y~x,family=binomial,data=a)) (M2=glm(y~x,data=b,family=binomial,weights=w)) However, with the arm and bayesglm function, the equivalent models with non-informative prior, give different results: library(arm) (M3=bayesglm(y~x,data=a,family=binomial,prior.scale=Inf, prior.df=Inf) )#M3=M2=M1 (M4=bayesglm(y~x,data=b,family=binomial,weights=w,prior.scale=Inf, prior.df=Inf))#M4 is different When I try a formulation with y=response rate or y=(response,failure), I get the rights coefs, but lower standard errors. c<-data.frame(y=c(4/9,6/7),x=c(4,5),w=c(9,7)) (M5=bayesglm(y~x,data=c,family=binomial,weights=w,prior.scale=Inf, prior.df=Inf)) y<-matrix(c(4,6,5,1),ncol=2) x<-matrix(c(4,5),ncol=1) (M6=bayesglm(y~x,family=binomial,prior.scale=Inf, prior.df=Inf)) Am I missing something? Thanks, Edouard Chatignoux, Statistician Health oservatory of the Paris Île-de-France region 43 rue Beaubourg- 75003 PARIS Tel. : 01 77 49 78 54 Fax. : 01 77 49 78 61 I'm running R 2.13.0 under windows XP > sessionInfo() R version 2.13.0 (2011-04-13) Platform: i386-pc-mingw32/i386 (32-bit) locale: [1] LC_COLLATE=French_France.1252 LC_CTYPE=French_France.1252 [3] LC_MONETARY=French_France.1252 LC_NUMERIC=C [5] LC_TIME=French_France.1252 attached base packages: [1] splines stats graphics grDevices datasets grid utils [8] methods base other attached packages: [1] foreign_0.8-44 arm_1.4-11 abind_1.3-0 R2WinBUGS_2.1-18 [5] coda_0.14-4 lme4_0.999375-39 Matrix_0.999375-50 lattice_0.19-26 [9] MASS_7.3-13 xtable_1.5-6 mgcv_1.7-6 ggplot2_0.8.9 [13] proto_0.3-9.2 reshape_0.8.4 plyr_1.5.2 loaded via a namespace (and not attached): [1] nlme_3.1-101 stats4_2.13.0 tools_2.13.0 ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
