I read Joseph Coveney's note in the Statalist of 2004 showing how to get
gllamm to estimate cross-classified models and this encouraged me to see
if I could get Stata to estimate a simple multiple membership model for
the purpose of spatial regression modelling of cancer statistics.
With a multiple membership model at the highest level only one variance
has to be estimated, but the random effect for a given observation is
weighted by the different memberships.
Unfortunately I've encountered a problem. The test data set I have
refers to observed and expected numbers of cancers at ward level all
located in a region, and the membership function is related to the
distance between each ward and the others. In this instance the
memberships all sum to unity but this needn't be so. MLWin has fitted a
Poisson model without covariates perfectly happily to these data with
the MCMC estimation option, and I want to compare the speed and output
of gllamm (while this may be slow, it may not be that bad as only one
3rd level variance is needed, and there's the convenience of a familiar
computing environment to consider).
The data & commands (adapted from Coveney's cross-classified model) are
given below. I'd naively assumed that the 3rd level intercepts could be
non-integer and represent the membership weights. However if run, the
result is a "too few variables specified" error message. Can anyone
spot what's going wrong?
clear
set more off
input ward observed expected distance lne w1 w2
w3 w4 w5 w6 w7 w8 w9 w10 w11
w12 w13 w14
1 5 12.473 8.645 2.524 .173 .104 .065
.074 .052 .066 .038 0 .058 .048 .094 .08
.09 .058
2 12 12.474 6.147 2.524 .078 .131 .099
.099 .084 .084 .062 .021 .068 .047 .085 .06
.047 .035
3 18 34.137 7.697 3.53 .044 .09 .119
.104 .102 .085 .045 .085 .062 .087 .065 .044
.041 .027
4 14 27.335 5.434 3.308 .05 .09 .105
.119 .1 .093 .051 .085 .061 .077 .058 .035
.035 .041
5 22 20.386 7.7 3.015 .033 .07 .093
.096 .109 .091 .056 .093 .072 .086 .071 .048
.049 .033
6 53 31.944 5.126 3.464 .042 .07 .082
.094 .092 .109 .071 .093 .075 .072 .064 .039
.047 .05
7 29 26.053 8.488 3.26 .027 .058 .077
.087 .095 .101 .122 .078 .087 .055 .063 .039
.059 .052
8 24 23.432 13.642 3.154 0 .022 .04
.052 .059 .071 .089 .137 .109 .108 .1 .07
.077 .066
9 26 26.422 6.061 3.274 .034 .054 .063
.074 .074 .089 .088 .067 .104 .073 .083 .062
.076 .059
10 21 19.028 8.36 2.946 .031 .041 .047
.06 .059 .075 .078 .082 .092 .115 .093 .074
.069 .084
11 11 29.72 1.811 3.392 .055 .066
.063 .075 .066 .081 .068 .047 .081 .072 .102
.08 .082 .062
12 28 19.802 5.701 2.986 .051 .051 .047
.061 .054 .071 .064 .057 .081 .089 .089 .109
.088 .088
13 36 37.129 8.488 3.614 .072 .049 .036
.052 .04 .061 .05 .044 .07 .083 .089 .107
.138 .109
14 30 28.049 10.072 3.334 .046 .035 .031
.047 .041 .062 .058 .066 .077 .1 .083 .109
.108 .137
end
/* NB w1-w14 denote the relative membership of wards 1-14 with each
other - a ward does not have membership of 1.00 with itself as the
memberships sum to unity */
set more off
/* Define the constant for the random intercept at level 2 */
generate byte ward_intercept = 1
eq k: ward_intercept
/* Assign each variable holding the multiple membership weights (to be
assigned to level 3 of the hierarchical model) to
an equation for a random intercept. NB the sum of w1..w14 = 1 */
eq r1: w1
eq r2: w2
eq r3: w3
eq r4: w4
eq r5: w5
eq r6: w6
eq r7: w7
eq r8: w8
eq r9: w9
eq r10: w10
eq r11: w11
eq r12: w12
eq r13: w13
eq r14: w14
/* Create a constant to identify level 3 */
generate byte region = 1
/* Define constraints so that a single variance will be estimated for
level 3 */
constraint define 1 w1 = w2
constraint define 2 w2 = w3
constraint define 3 w3 = w4
constraint define 4 w4 = w5
constraint define 5 w5 = w6
constraint define 6 w6 = w7
constraint define 7 w7 = w8
constraint define 8 w8 = w9
constraint define 9 w9 = w10
constraint define 10 w10 = w11
constraint define 11 w11 = w12
constraint define 12 w12 = w13
constraint define 13 w13 = w14
/* Estimate the model */
gllamm observed distance, i(ward region) family(poisson) offset(lne)
nrf(1, 14) eqs(k r1 r2 r3 r4 r5 r6 r7 r8 r9 r10 r11 r12 r13 r14)
nocorrel ///
constraints(1 2 3 4 5 6 7 8 9 10 11 12 13) nip(3) adapt trace allc
end
Paul Silcocks BM BCh, MSc, FRCPath, FFPH, CStat Clinical Senior Lecturer
Nottingham Clinical Trials Support Unit
Room B39
School of Community Health Sciences
University of Nottingham Medical School
Nottingham
NG7 2UH
tel +44 115 9515151 x 30514
Direct line: +44 115 8230505.
General office: +44 115 8230500
fax +44 115 8230515
e-mail paul.silcocks@nottingham.ac.uk
This message has been checked for viruses but the contents of an attachment
may still contain software viruses, which could damage your computer system:
you are advised to perform your own checks. Email communications with the
University of Nottingham may be monitored as permitted by UK legislation.
*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/