Hello All
Here is a program which models injury rates with spatial variation in
WinBugs.
Here is a summary of the model,
b[i] represent the CAR distribution for spatial effects in the model
Weights are assigned with a for loop for the CAR.normal distribution
expl[i] is an 83x1 matrix and the product of beta * X, where beta is a 11
by 1 vector of unknown parameter coefficients and X is a 83 by 11 matrix of
11 explanatory variables for 83 local health areas(LHA).
n[i] has a Poisson distribution and n[i] contains the number of injuries in
each LHA/100,000 persons.
log(lambda[j]) is the expected rate with population(p[j] in LHA as offset
and insludes the expl[j] effects, b[j] spatial effects, alpha0 is a
convolution prior for CAR.normal statement, and alpha1 a random effects
variable to account for overdispersion.
the prior distributions are chosen to in non-informative
----------------------------------------------------------------------------
I am concerned with my mathematical reasoning for this model. I want to map
the injury rates for each LHA. I assume that this value is lambda[i].
Please comment on whether I obtain the correct injury rates, if at all,
from this code. The entire code for the model is below. The model compiles,
updates, and runs without any errors.
I did notice that the alpha0 and alpha1 density graphs are all over the
place as well as the hisotry for lambda[83] does not seem right.
Please contact me ASAP with any comments.
Thanks
Ayaz
#MODEL.
model;
{
b[1:N] ~ car.normal(adj[],weights[],num[],tau)
for(k in 1:sumNumNeigh){
weights[k]<- 1
}
for (i in 1:N){
expl[i]<- beta1*X[i,1]+beta2*X[i,2]+beta3*X[i,3]+beta4*X[i,4]
+beta5*X[i,5]+beta6*X[i,6]+beta7*X[i,7]+beta8*X[i,8]+beta9*X[i,9]+beta10*X
[i,10]+beta11*X[i,11];
}
for (j in 1:N) {
n[j] ~ dpois(lambda[j])
log(lambda[j]) <-log(p[j])+expl[j]+b[j]+alpha0+alpha1;
}
# non-informative priors.
alpha0 ~ dflat()
tau ~ dgamma(0.5,0.0005)
alpha1 ~ dnorm(0,0.00001)
beta1~ dnorm(0,0.000001)
beta2~ dnorm(0,0.000001)
beta3~ dnorm(0,0.000001)
beta4~ dnorm(0,0.000001)
beta5~ dnorm(0,0.000001)
beta6~ dnorm(0,0.000001)
beta7~ dnorm(0,0.000001)
beta8~ dnorm(0,0.000001)
beta9~ dnorm(0,0.000001)
beta10~ dnorm(0,0.000001)
beta11~ dnorm(0,0.000001)
sigma <- sqrt(1/tau);
}
----------------------------
# Model data.
list(N=83,sumNumNeigh=370,
num=c(
2,4,4,5,4,7,5,5,6,3,
5,5,3,6,4,5,4,8,7,3,
12,6,6,4,6,10,4,4,8,7,
6,4,4,3,4,1,0,1,2,3,
7,6,3,2,2,3,11,4,0,7,
4,5,5,7,8,8,2,5,2,4,
2,0,3,5,3,4,3,6,3,8,
4,6,4,4,5,2,0,3,4,5,
6,3,2),
adj=c(
2,4,
1,3,4,5,
2,4,5,6,
1,2,3,6,17,
2,3,6,7,
3,4,5,7,9,17,18,
5,6,8,9,10,
7,9,10,11,21,
6,7,8,11,18,21,
7,8,11,
8,9,10,12,21,
11,13,14,21,22,
12,14,15,
12,13,15,21,22,73,
13,14,16,73,
15,22,30,31,73,
4,6,18,56,
6,9,17,19,21,25,56,74,
18,20,21,23,25,30,74,
19,21,74,
8,9,11,12,14,18,19,20,22,23,30,74,
12,14,16,21,30,73,
19,21,24,25,29,30,
23,25,26,29,
18,19,23,24,26,56,
24,25,27,28,29,48,54,56,70,79,
26,54,55,56,
26,29,47,70,
23,24,26,28,30,31,47,72,
16,19,21,22,23,29,31,
16,29,30,32,47,72,
31,33,71,72,
32,34,41,71,
33,35,41,
34,36,41,42,
35,
40,
40,42,
38,39,42,
33,34,35,42,47,71,72,
35,39,40,41,43,47,
42,44,47,
43,47,
46,47,
45,47,70,
28,29,31,41,42,43,44,45,46,70,72,
26,54,75,79,
51,52,55,80,81,83,82,
50,75,81,82,
50,53,54,55,81,
52,54,55,75,81,
27,26,48,52,53,55,75,
27,50,52,54,53,56,58,80,
17,18,25,26,27,55,57,58,
56,58,
55,56,57,76,80,
60,61,
59,61,63,64,
59,60,
60,64,65,
60,63,65,66,68,
63,64,66,
64,65,67,68,
66,68,69,
64,66,67,69,70,78,
67,68,70,
26,28,46,47,68,69,78,79,
32,33,41,72,
29,31,32,41,47,71,
14,15,16,22,
18,19,20,21,
48,51,53,54,81,
58,80,
68,70,79,
26,48,70,78,
50,55,58,76,83,
50,51,52,53,75,82,
50,51,81,
50,80),
n=c(0, 1,0, 1, 0, 0, 2, 0, 0, 1, 1, 3, 3, 1, 1, 1, 1, 0, 4, 1, 6, 5, 8, 2,
2, 6, 4, 0, 1, 0, 2, 5, 8, 10, 20, 3, 6, 12, 2, 2, 6, 4, 3, 3, 2, 0, 4, 0,
2, 0, 0, 0, 2, 3, 4, 14, 3, 4, 6, 3, 4, 0, 5, 0, 0, 7, 0, 6, 7, 3, 1, 2, 0,
0, 2, 0, 2, 0, 5, 2, 4, 0, 0),
p=c(20857,31224, 9041, 9848, 13933, 3810, 28129, 15244, 5545, 24989, 9081,
4279, 15606, 38421, 4468, 5273, 10221, 11390, 31887, 11397, 62950, 144817,
114602, 18196, 7140, 38953, 36020, 6931, 9353, 15844, 10410, 83150, 136190,
134943, 393605, 120576, 157757, 852802, 38597, 164424, 94422, 207350,
131805, 41641, 28660, 22941, 30835, 5030, 7635, 1434, 28498, 9499, 26585,
11590, 26766, 142078, 42791, 42952, 177254, 65893, 55909, 12105, 60744,
7845, 17636, 98859, 32394, 41568, 65535, 49914, 48539, 8384, 11290, 8343,
18903, 9707, 2949, 5554, 23305, 2198, 33725, 3356, 1016),
X=structure(.Data=c(
0.61171001,-0.34845999,-1.4647501,-0.38999999,-0.49979001,-0.44328001,-
0.33465999,0.53478998,-0.23999999,-0.68835002,-0.04656,
0.36151001,-0.47507,0.95192999,-0.14,-0.31869,-0.16855,-
0.24770001,0.30066001,1,-0.39714,-0.24216001,
0.31341001,-0.51988,-0.28128999,0.47,0.2287,-0.51587999,-0.47696,-0.28033,-
0.36000001,-1.0626,-1.20188,
0.73925,-0.51345998,-1.85254,-0.33000001,-0.10374,-0.53551,-
0.55075002,1.54918,0.09,0.71414,1.4895101,
1.22129,-0.48984,-1.32783,1.45,1.91557,-0.54531997,-0.58236998,2.0387199,-
0.23,-0.95172,0.46193999,
0.03657,-0.32482001,-1.7947,0.94,0.64238,-0.56298,-0.64824998,1.7243201,-
1.35,-0.97795999,-1.20188,
-0.56715,-0.35563001,1.79805,-0.25,0.17296,-0.52372998,-0.41372001,-
0.62015998,1.15,-0.66961998,-0.52869999,
0.61862999,-0.54323,-0.59566998,0.03,0.78408003,-0.51587999,-0.47433001,-
0.23710001,0.40000001,-0.30131999,-0.44784001,
0.15538999,-0.01742,-0.37066999,1.39,1.41031,-0.52961999,-
0.50331002,0.21235,-0.58999997,-0.45190999,-0.87538999,
0.64012003,-0.43952,-0.79883999,-0.2,0.75149,-0.29806,-
0.23980001,0.25995001,3.5599999,-0.29438999,0.27351001,
0.5952,-0.50357997,0.11315,0.97000003,0.78596997,-0.55120999,-
0.47960001,1.9220901,0.05,0.93150002,0.25457001,
0.18925001,-0.37537,1.03028,2.0699999,-1.03101,-0.54923999,-
0.57182997,0.94300002,-2.1800001,-0.72612,-0.23813,
1.54142,0.06205,-1.2282701,1.75,1.3811001,-0.35497001,-0.49541,-0.44477999,-
1.29,-1.0024,-0.86149001,
0.26464999,-0.01385,1.17828,0.94999999,0.56357002,-0.24704,-0.13966,-
0.86594999,2.0999999,0.19354001,-0.27158001,
0.93395001,-0.14215,0.12647,2.49,0.73056,-0.56493998,-0.68777001,-
0.40656999,-1.99,-0.07531,0.22733,
0.55800003,-0.40838999,-1.08433,1.1900001,2.57356,-0.55905998,-
0.59290999,0.26563001,-0.44,0.73940998,-0.21128,
0.01686,-0.31450999,-0.32956001,0.12,-0.66896999,-0.41973001,-
0.55338001,2.63538,-0.15000001,0.42096001,-0.47220001,
0.16143,-0.41604,-0.06375,0.38999999,-0.92866999,-0.51392001,-
0.44007,0.76475,-0.33000001,0.43507999,0.58397001,
1.14299,-0.32376999,-0.46697,0.94999999,-0.09758,-0.42561999,0.01054,-
0.21979,0.02,-0.56696999,-0.98540002,
1.0707901,-0.44431001,-0.17399,0.22,-0.78807002,-0.54139,-
0.57446003,1.3516001,-1.04,-1.17344,-1.0218199,
0.30023,-0.31819999,1.37573,0.54000002,-0.30001,1.72904,0.31358001,-
0.02333,1.37,-0.68114001,-0.71648002,
0.43299001,-0.29604,0.65854001,-0.11,-
0.54236001,0.55948001,1.66278,0.57283002,1.0599999,-0.39825001,-0.30807999,
-0.09898,-0.27928999,1.35982,-0.23,-
0.35743999,0.52023,0.92229998,0.21343,0.88,-0.35732001,-0.58174998,
0.74400997,-0.12174,-1.44078,1.6900001,0.90358001,-0.48449001,-
0.59026998,0.03344,-1.0700001,-0.13698,0.40742999,
0.24898,-0.35047001,-2.1259501,0.88999999,0.32874,-0.57082999,-0.68777001,-
0.42085001,-1.64,0.33349001,0.25997999,
-0.06203,-0.35745001,-0.75247997,-0.12,0.68519998,-0.31376001,-
0.47169,0.46496999,-0.25,0.95262998,0.57204998,
0.03558,-0.30371001,0.31595001,0.75,0.23944999,-0.54728001,-
0.45852,1.74463,-0.55000001,0.39456999,0.11656,
-2.2797201,-0.46574,1.4936399,0.80000001,1.67726,-0.56298,-
0.57972997,1.66293,-0.44999999,1.22694,1.1055599,
-0.71847999,-0.24150001,-0.06511,1.27,0.41427001,0.10029,-
0.42688999,0.13407999,-1.34,2.5676799,2.5945201,
-0.45062,-0.33227,1.01995,1.23,1.66254,-0.48644999,-0.57446003,0.46158999,-
1.25,1.44692,1.13193,
-0.26721999,-0.18332,1.10382,1.6900001,0.53241998,-0.53157997,-
0.45852,0.14308999,-1.09,2.4198501,0.62128001,
0.56379998,-0.27809,0.89507997,0.40000001,-0.50215,-0.38047999,0.71938998,-
0.63069999,0.41,-0.70354998,-0.60677999,
1.09961,0.75323999,-0.19119,-0.34999999,0.05278,0.73804998,1.5573699,-
0.01684,0.05,-0.50340003,-0.41668999,
0.86539,-0.07916,-0.42824,-1.17,-0.90676999,0.54574001,1.90785,-0.26817,-
0.25999999,-0.96853,-0.84960002,
0.59601003,1.67006,0.49090001,-0.97000003,-0.64363003,2.2471001,4.7485299,-
0.60598999,-0.25,-0.82212001,-0.80101001,
1.15658,0.68307,-0.83095998,-2.5,-1.21052,0.30634001,0.96183002,-
0.88565999,0.36000001,-0.01092,-0.68769002,
0.47426999,5.5380001,-0.36542001,-2.3399999,-1.48518,-
0.12537999,0.86168998,-0.49015,0.76999998,-1.00457,-0.80954999,
-2.2369399,3.3706601,0.59042001,0.67000002,-1.16993,1.9959199,4.8196802,-
0.55467999,6.3099999,-0.94669002,-0.88024998,
-2.1933,1.33276,1.29317,-1.12,0.09984,-0.42954001,-0.11068,-
0.36348999,8.0600004,-0.26082,-0.5643,
-0.83041,3.8418801,0.48895001,-1.67,-1.08471,2.5257499,1.5573699,-
0.62199003,0.56,-1.13423,-0.81010997,
0.66140002,0.18894,0.00989,-0.80000001,-0.50120997,-0.22153001,0.66933,-
0.21841,-0.15000001,-0.58165997,-0.61813998,
0.4129,2.10832,0.44538,-2.26,-1.42214,1.26004,1.48622,-0.11695,-0.34,-
0.58702999,-0.60035998,
-0.14715999,1.49322,0.16819,-3.8,-1.31708,0.40838,0.64824998,-
0.49250001,1.76,-1.21551,-0.90551001,
0.74619001,1.86726,-1.21259,-6.1900001,-1.22936,0.26905,0.07378,-
0.88647002,2.6199999,-0.28053999,-0.47492999,
0.39443001,-0.37086001,0.5291,-0.30000001,-0.30577001,-0.39225999,0.01318,-
1.06761,0,-0.51090002,-0.10985,
0.42860001,-0.45436999,-0.0843,0.40000001,0.29697001,-0.46094,-0.40318,-
0.64411002,0.66000003,-0.47497001,-0.45284,
-1.35882,0.27480999,-0.48291001,-2.3599999,-0.92860001,-0.25685,-
0.20291001,-1.04248,-0.69999999,-0.01982,-0.12761,
-0.98080999,-0.47801,-0.31740999,0.04,1.17063,-0.57082999,-0.68777001,-
1.26564,0.27000001,0.42157,-0.59241998,
-1.29354,-0.06236,-0.82713997,-0.01,-1.03238,-0.34711999,-0.59290999,-
0.66509998,-0.36000001,1.41466,2.81442,
-0.12845001,-0.30952001,-1.57569,0.63,-1.52492,-0.56690001,-0.61926001,-
0.54234999,2.6099999,2.3141,0.58064997,
-1.16315,-0.09961,1.26999,-0.04,0.63665998,-0.03511,-0.46114999,-
0.30842,0.20999999,1.04453,0.75538999,
-1.99599,-0.43557999,0.77165997,1.1900001,1.01816,-0.50410998,-0.63770998,-
0.52181,-0.86000001,1.91346,-0.10985,
0.19772001,-0.37103,-1.43538,-0.79000002,-0.34371999,-0.40206999,-
0.31621999,-0.40493,0.22,0.11037,0.43981001,
0.18942,-0.29837999,-0.39857,0.02,1.42478,-0.55317003,-
0.62980002,1.4873199,-0.72000003,0.05252,-0.69700998,
-0.09725,-0.40490001,-0.73532999,0.62,0.32723999,-0.49822,-0.41635001,-
0.74892998,-1.3,1.10605,2.6784999,
-0.33153,-0.33101001,0.79668999,-0.70999998,-
0.69098002,2.10973,1.21217,0.10119,0.5,0.45016,-0.22269,
-0.03411,-0.43299001,0.6803,-0.16,-0.59798002,-0.23723,-0.01054,2.1036999,-
0.2,1.57481,0.84815001,
-0.06519,-0.27278,-0.73895001,-0.37,-1.1246001,-0.51392001,-
0.08959,1.54697,-1.36,-0.22576,0.00162,
-1.54645,0.12261,1.5632499,-1.79,-0.29513001,1.47001,0.47433001,-
0.50739002,5.3699999,-0.82810003,-0.55954999,
0.80699998,-0.21162,0.48337999,-0.95999998,-1.00229,0.7773,0.06324,-
1.07858,-3.0999999,-0.33706,-0.55958003,
1.76157,-0.2554,-1.64019,-2.5799999,-0.31511,2.7553501,1.46514,-0.44881999,-
0.81999999,0.81114,0.004,
0.29943001,0.03257,-0.54157001,-1.8099999,-0.23443,-0.57082999,-0.68777001,-
0.81090999,0.41,-1.14261,-0.98038,
0.70750999,-0.42822,0.04891,-0.1,-0.64937001,0.63600999,0.05007,-
0.44782001,0.12,-0.61229002,-0.30226001,
0.00028,-0.35982001,1.7847,1.33,-1.22527,-0.36085999,-0.60872,-0.26063001,-
1.8,-0.60861999,-0.91658002,
0.73219001,-0.36353999,-0.42625001,0.30000001,1.0229501,-0.48449001,-
0.55075002,-0.21051,-0.49000001,-0.1218,0.02247,
-0.23308,-0.22394,1.88796,0.01,-0.33324,1.44842,0.75102001,-0.14732,1.49,-
0.40166,-0.35416999,
1.3821,-0.34740999,-0.70656002,0.51999998,-0.34551001,-0.37263,-0.31885001,-
1.04882,-0.38999999,0.25137001,-0.53702003,
-0.19253001,-0.4149,0.99971998,0.74000001,-0.242,-0.45504999,-0.10804,-
1.21986,-0.63999999,0.20468999,0.77569002,
0.42642999,-0.3238,0.82628,0.12,-0.28246999,-0.40403,-0.10277,-
0.35049,0.95999998,-0.62311,-0.48899999,
-0.07141,-0.35132,0.88905001,-0.16,-0.58157998,-0.34908,-
0.07905,1.2898,0.1,-0.30223,-0.14122,
0.46419999,-0.17621,0.32374001,0.38,-0.36182001,-0.16855,-0.0527,-1.00927,-
0.91000003,-0.77749997,-0.89573997,
-0.18038,-0.06386,0.81682998,0.91000003,-0.97605002,-0.45113,-0.41635001,-
1.18189,-2.73,-0.75660998,-0.40858001,
1.66392,-0.34915,-0.77809,0.51999998,0.66942,-0.37852001,-0.25034001,-
0.90526003,0.06,-1.01565,-0.23792,
0.46912,-0.58117002,1.2506599,2.1500001,1.1024899,-0.53355002,-
0.15019999,0.91703999,-0.54000002,0.54118001,1.16298,
0.44869,-0.36067,-0.09387,-0.66000003,0.35354999,5.4418001,-0.49013999,-
1.25164,-0.19,-0.61786002,-0.46785,
-1.06415,-0.32424,-1.62603,-0.67000002,0.56888998,-0.48644999,-
0.50858003,2.8912101,-0.72000003,2.10935,2.4042799,
-2.0711901,-0.36937001,0.36980999,2.46,3.9033899,-0.57082999,-0.68777001,-
1.36797,-0.92000002,-1.58085,-1.20188,
0.77574003,-0.45104,-2.0343599,-1.27,-1.4738801,-0.41973001,-0.61926001,-
1.00287,-1.48,-0.19595,0.55295002,
-0.73705,-0.25294,-0.02287,-0.61000001,0.85367,-0.33730999,-
0.46641999,0.09577,-1.01,0.65726,1.3650301,
-2.6903,-0.62396002,1.2036999,-0.75999999,-0.86066997,-0.57082999,-
0.68777001,-0.472,-1.3200001,2.61834,3.0929501,
-0.41396999,-0.36122,0.38883999,-0.56,-0.56156999,-0.42954001,-0.35574001,-
0.48331001,0.28999999,0.35251999,0.18927,
-0.84397,-0.62396002,-0.79931998,2.5,0.82523,-0.57082999,-
0.68777001,2.08915,-2.6400001,2.1514001,3.01512,
-3.3134,-0.62396002,0.71851999,1.46,0.68796003,-0.57082999,-0.68777001,-
1.36797,-1.84,-1.7903301,-1.20188),.Dim=c(83,11)
)
------------------
#Initial data.
list(alpha0=0,beta1=0, beta2=0, beta3=0,beta4=0,beta5=0,beta6=0, beta7=0,
beta8=0, beta9=0, beta10=0, beta11=0,
alpha1=0,tau=1,
b=c(
0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0))
)
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
