Saturday 4 August 2007

Archaeologists Beware!


Here we have the latest technology - a trap to catch archaeologists:



How does it work, you are asking yourselves. Well, it's based on the observation that archaeologists are fascinated by old rubbish. So, the trap is set up to mimic a mediaeval rubbish bin. The archaeologist is immediately attracted, and leaps in. When he is in there up to the middens, the trapper leaps out and ties up the archaeologist.

The one question I cannot answer is why anyone would want to capture archaeologists. Pest control?
Read more!

Friday 3 August 2007

stRing

Here is some more code, to keep the troll confused and happy.


library(ggm)

# Figure 3.1 data
Data=read.table("http://pages.usherbrooke.ca/jshipley/recherche/NINA%20files/fig31.dat")

# Write the DAG
DAG1=DAG(V2~V1, V3~V2, V4~V2, V5~V3, V5~V4, V5~V6)
# Draw the DAG
drawGraph(DAG1)
# D-separation test
dSep(DAG1, "V3","V1", NULL)
# Shipley Test
shipley.test(DAG1, cov(Data), dim(Data)[1])

#####################################################
#### SEM

Data=read.table("http://pages.usherbrooke.ca/jshipley/recherche/NINA%20files/fig31.dat")
names(Data)=LETTERS[1:dim(Data)[2]]

library(sem)
# Input the model
path3.1=specify.model()
A -> B, beta.AB, NA
B -> C, beta.BC, NA
B -> D, beta.BD, NA
C -> E, beta.CE, NA
D -> E, beta.DE, NA
F -> E, beta.FE, NA
A <-> A, var.A, NA
B <-> B, var.B, NA
C <-> C, var.C, NA
D <-> D, var.D, NA
E <-> E, var.E, NA
F <-> F, var.F, NA

# Fit the SEM
model3.1=sem(path3.1, cov(Data), N=dim(Data)[1])
summary(model3.1)

# Specify an alternative model: remove the B->D path
path3.2=specify.model()
A -> B, beta.AB, NA
B -> C, beta.BC, NA
C -> E, beta.CE, NA
D -> E, beta.DE, NA
F -> E, beta.FE, NA
A <-> A, var.A, NA
B <-> B, var.B, NA
C <-> C, var.C, NA
D <-> D, var.D, NA
E <-> E, var.E, NA
F <-> F, var.F, NA

path3.2

model3.2=sem(path3.2, cov(Data), N=dim(Data)[1])
summary(model3.2)

# Compare the two models
anova(model3.1, model3.2)

###################
# String

# Read data (use your own file!)
string=read.table("string.txt")

# First parameterisation
path.string1=specify.model()
Length -> V1, NA, 1
Length -> V2, beta.2, NA
Length -> V3, beta.3, NA
Length -> V4, beta.4, NA
V1 <-> V1, var.1, NA
V2 <-> V2, var.2, NA
V3 <-> V3, var.3, NA
V4 <-> V4, var.4, NA
Length <-> Length, var.L, NA

model.string1=sem(path.string1, cov(string), N=dim(string)[1])
summary(model.string1)

# First parameterisation, Var(V1) set to zero
path.string1a=specify.model()
Length -> V1, NA, 1
Length -> V2, beta.2, NA
Length -> V3, beta.3, NA
Length -> V4, beta.4, NA
V1 <-> V1, NA, 0
V2 <-> V2, var.2, NA
V3 <-> V3, var.3, NA
V4 <-> V4, var.4, NA
Length <-> Length, var.L, NA

model.string1a=sem(path.string1a, cov(string), N=dim(string)[1])
summary(model.string1a)

# Second parameterisation
path.string2=specify.model()
Length -> V1, beta.1, NA
Length -> V2, beta.2, NA
Length -> V3, beta.3, NA
Length -> V4, beta.4, NA
V1 <-> V1, var.1, NA
V2 <-> V2, var.2, NA
V3 <-> V3, var.3, NA
V4 <-> V4, var.4, NA
Length <-> Length, NA, 1

model.string2=sem(path.string2, cov(string), N=dim(string)[1])
summary(model.string2)

# Second parameterisation, Var(V1) set to zero
path.string2a=specify.model()
Length -> V1, beta.1, NA
Length -> V2, beta.2, NA
Length -> V3, beta.3, NA
Length -> V4, beta.4, NA
V1 <-> V1, NA, 0
V2 <-> V2, var.2, NA
V3 <-> V3, var.3, NA
V4 <-> V4, var.4, NA
Length <-> Length, NA, 1

model.string2a=sem(path.string2a, cov(string), N=dim(string)[1])
summary(model.string2a)

# Profile likelihood, to look at likelihood of Var(V1)
GetDev=function(var, path, Cov, N, ...) {
path[5,3] = var
deviance(sem(path, S=Cov, N=N, ...))
}

Var1=as.character(seq(0,5, length=100))
Deviance=sapply(Var1, GetDev, path=path.string1a, Cov=cov(string), N=dim(string)[1])

plot(Var1, -Deviance+min(Deviance), type="l")
abline(h=-3.71, col=2)

DevDiff=-Deviance+min(Deviance)
max(Var1[which(DevDiff>-3.71)])
Read more!

Thursday 2 August 2007

How long is a piece of string?


Possibly the best post I will ever make....

43.08322 5.633989 14.39994 43.0198
37.27808 -0.6305197 13.60435 35.97464
56.98669 5.910917 20.24388 65.7003
37.9433 -6.700168 16.72924 33.24414
41.87214 4.740604 18.91497 21.96089
55.57904 -7.990423 27.50937 59.05561
67.64643 11.91912 29.99546 57.90356
44.71473 -5.295052 12.39872 43.90684
37.87198 5.177735 16.46788 36.11007
51.48916 5.99543 22.11037 46.40012
39.26152 15.24808 14.50037 33.99699
56.37837 -14.86887 18.16238 54.52073
39.0634 1.106713 15.83562 35.57628
59.4622 15.60182 22.19649 53.731
51.79565 -2.257035 21.14818 40.94382
49.0601 1.274861 16.57190 60.04407
45.4396 1.187506 23.36360 32.59971
24.84383 0.02483562 19.02628 27.39255
42.39041 10.86633 16.92249 31.95631
38.49579 3.773862 8.178377 35.90047
54.93457 6.60548 15.48205 41.43354
53.47267 12.77436 21.74509 68.12006
45.73146 19.44002 19.95146 41.04718
30.45508 15.49729 9.07437 36.12982
31.09868 15.89208 14.05893 25.72002
55.53552 18.22403 23.28495 49.00676
53.61962 8.500991 22.92718 43.00829
57.66102 -5.72349 26.62879 54.14032
44.29104 -5.413569 14.56123 29.56197
27.98142 12.27632 7.381431 33.86140
39.29747 5.73134 13.61576 32.31704
63.61147 -4.124622 20.77645 65.57825
48.30483 -0.9741744 23.14768 64.53353
39.40313 1.301428 8.893243 34.26687
51.05807 7.808775 22.74802 56.97688
54.37367 1.224387 16.43467 44.20504
37.14993 7.237878 12.33427 42.1512
34.96905 -9.276209 13.83203 34.83956
64.49703 -4.036688 30.76852 59.17002
57.57922 4.81336 19.08255 58.9137
42.89432 -3.014686 15.16899 51.37077
50.47695 1.228656 20.90042 52.62307
60.55926 12.58458 19.93596 50.11331
44.0029 10.0013 18.27015 45.40107
42.21744 6.606347 22.75961 58.13801
50.91977 18.42849 21.01416 49.05971
43.72696 10.20114 16.31195 57.22922
60.6012 1.190260 26.42322 70.65941
65.08044 11.42588 23.49145 63.92476
68.47565 -1.136795 25.78938 87.5605
46.33998 7.442789 20.58232 59.85457
50.79214 6.526056 19.52277 48.91348
45.64809 3.005203 17.14066 43.24646
41.45269 5.170724 11.35216 36.99929
53.06 11.57520 18.53381 54.74449
56.54957 9.515325 26.00357 58.91207
57.34748 4.723733 18.30328 54.43463
50.26852 9.466306 21.77719 53.3892
31.54898 -4.326571 8.553822 35.21168
62.9613 8.670223 22.37621 48.70231
51.32045 2.452138 21.85014 44.80914
48.21747 7.79318 25.90706 47.41505
48.82583 14.84098 13.47211 53.04329
57.1413 11.25805 18.08639 49.65892
53.25459 -4.27024 18.02547 68.7887
55.97105 6.01462 26.85532 61.14111
46.90504 14.30717 15.42983 56.82181
67.18347 0.3422343 28.47143 53.90871
60.51425 18.29886 23.32582 59.69433
58.43662 6.812244 18.30054 50.66552
55.34559 4.416507 20.75887 56.4903
51.46907 6.437716 23.76091 46.35918
55.11886 7.143447 18.90608 61.75727
63.13792 9.285244 20.93044 61.66314
61.87049 7.384118 26.37211 40.59172
74.20589 7.355128 32.55491 87.44812
55.10594 -0.1525557 21.72056 51.75426
24.05928 5.307248 7.631405 14.88209
51.80546 4.658656 20.88394 47.48967
37.46283 -10.4277 11.45324 21.98236
40.32689 7.986716 13.30132 54.67371
38.10324 1.991273 12.84804 55.30039
55.2736 10.60343 18.49757 38.21288
47.07817 20.04340 16.51424 49.3371
58.96848 5.799431 20.01643 55.83296
71.49642 9.116315 27.19482 70.67184
62.60924 -1.557933 28.71131 70.46987
51.13872 11.74212 17.35272 36.62495
39.0086 2.329680 19.28133 34.90509
76.80814 16.32340 29.29323 77.37185
44.05229 18.07265 13.12489 37.31700
44.58876 17.50309 17.30197 40.72883
66.08837 9.521458 24.62508 67.92723
44.70895 -9.616063 24.35619 40.39696
48.95169 11.08539 16.25009 46.68092
72.20372 6.56739 26.02335 65.54035
36.72116 6.479657 12.10085 46.69272
44.6516 -0.4680244 19.14367 48.71209
46.10673 -8.10203 12.65460 45.67924
56.1415 11.15108 21.70855 74.11605
46.78581 9.729058 18.70586 54.91059
49.20233 -0.2790261 20.29539 54.23667
45.62776 15.40980 21.89740 36.38535
42.99238 6.919789 15.50135 47.28677
54.43287 12.92555 16.94803 52.39224
44.51279 1.02478 22.49600 52.49208
43.64437 1.412202 10.43212 53.15029
34.8613 8.567357 12.73069 39.31951
57.87542 0.9947726 22.30917 67.91817
55.27183 9.722042 18.51507 78.40459
31.71783 -7.593688 12.03276 43.22243
46.90415 11.22460 22.13732 48.0647
17.89495 -1.989388 8.589473 32.64595
47.26154 18.13389 20.93814 50.39569
51.97978 -0.742969 21.92611 61.26977
37.6244 1.130264 13.08287 42.63462
57.1765 -1.871521 14.04430 68.6356
40.156 -8.217474 15.43289 25.20931
66.52183 -0.3503514 27.62605 48.65772
44.29208 7.897813 25.55882 30.89449
64.3822 -1.209621 23.14808 54.64099
59.56115 9.57178 20.09612 65.67535
47.23148 2.225237 15.67632 51.25574
47.78437 5.957391 23.0804 31.84572
51.71551 14.18624 18.76672 42.09963
50.41613 -2.764840 18.76497 40.23106
57.31479 5.689861 26.50951 60.05792
41.27538 6.095382 12.22231 39.24507
34.07432 3.550837 15.19683 28.38449
51.88646 -0.2055126 21.60829 53.51454
61.70346 -21.06496 21.50972 63.80856
57.57247 -1.213927 23.48830 53.32069
49.67057 2.958423 14.83082 61.50704
54.14113 8.83217 23.47951 48.58643
49.89234 16.73157 18.73319 32.56289
53.17817 -1.315715 23.10059 59.69463
67.96876 18.46018 29.49493 70.76115
34.20969 8.591464 15.65683 36.51379
65.76896 7.56501 26.12036 61.40985
62.34774 14.99913 24.00971 41.68162
47.74665 3.605466 20.40071 17.66518
60.26479 1.957539 21.96315 50.42938
47.2237 4.228247 22.50998 35.00704
45.04607 -8.937676 16.75411 53.25114
42.14120 4.965197 15.69227 39.63621
42.07501 3.973461 18.43234 48.03502
52.67104 7.42322 16.60783 48.78771
45.87456 1.693674 25.34206 53.88313
38.80739 -1.168698 14.54299 37.1439
60.12847 10.67939 18.64974 44.68634
58.28223 6.430793 22.32087 46.62149
43.454 9.3688 18.34806 60.74563
61.15555 6.758348 24.55251 61.19742
55.80496 6.91809 22.22007 41.10921
50.49143 -0.1087601 19.30260 51.63405
57.61026 2.042491 14.69968 54.28266
41.55212 -0.1718027 17.34501 34.83176
59.52748 -0.6971561 23.90598 62.33222
48.58037 13.17525 22.31383 54.61252
36.16443 3.855484 10.41436 39.14067
34.8746 13.68047 9.348352 37.23884
45.38424 -9.679704 16.41429 50.2606
38.01672 -5.733999 14.78290 38.90471
50.74452 4.605116 19.33822 56.48359
36.71898 9.331181 8.60925 33.40870
52.15962 2.773299 24.30411 39.83019
73.61568 -0.7022454 23.0522 64.07831
48.57496 8.681073 18.71645 38.35804
35.09403 -3.317196 13.67876 35.55393
41.14362 0.2859828 16.03258 55.41728
60.42156 19.86611 25.99627 60.13635
57.91767 18.09433 23.74566 38.92939
68.5777 2.792227 27.17936 64.1104
38.67721 2.394332 14.74335 51.11878
43.33151 6.654543 18.17277 34.78806
32.16489 5.378586 17.69214 38.28201
38.09104 2.722748 17.04729 38.5717
57.65442 -1.166543 27.38046 63.76525
42.85036 -0.1664495 19.76169 42.86057
30.81336 -18.98881 17.1611 49.82694
44.81933 -0.1029979 14.00432 35.39499
49.04004 9.404016 23.09975 44.44873
33.7447 -1.342813 15.80708 10.53570
38.79723 3.765855 21.6306 24.81213
66.85159 19.75231 19.79586 53.1625
26.79097 0.817258 11.02373 56.88012
37.20261 0.173163 15.84341 23.58919
33.37625 -8.135147 8.814127 27.44279
35.96085 -1.496246 17.03312 16.22190
34.84232 -2.665887 12.31025 38.5161
46.4912 0.1855528 19.53513 37.71706
39.82818 -1.251567 15.27734 36.08796
42.90098 -11.78600 17.27517 51.89914
51.62899 10.44025 18.37046 57.14262
41.70951 13.13171 17.49616 37.1912
47.13227 -1.742686 23.88572 49.90174
32.02646 5.165141 10.78761 40.83208
58.07391 -1.254475 28.28493 61.74809
38.52684 -1.899376 11.82144 45.43454
64.9634 4.697399 18.85557 60.94828
42.33774 -5.140385 18.92792 33.66515
40.09141 5.735004 12.07397 30.23979
36.77622 -1.373729 8.959085 50.2355
48.24848 -2.477777 26.38349 34.84796
46.80242 3.878493 14.63697 53.24456
60.48241 1.329867 23.34244 60.8137
39.81684 2.536786 15.46595 45.34798
56.37617 9.731282 23.93499 60.20108
71.21766 1.948052 28.54894 76.46366
61.44952 0.03265491 24.71560 79.16348
59.33382 10.14008 20.35875 64.21955
49.90094 7.110765 22.86384 46.81132
57.85128 2.27204 21.54030 62.51414
56.35488 0.02495596 23.93359 50.76844
51.82123 1.897808 23.42504 47.90481
39.87752 3.817866 24.09952 36.94029
53.04694 0.943916 19.43468 66.65729
69.87825 2.468058 26.86654 80.04521
49.20798 8.266817 22.22331 44.30947
32.68928 3.734302 14.28460 26.69898
44.14841 -1.230833 17.73934 47.07987
46.62789 9.472355 22.57244 47.83675
58.20436 10.32852 20.23259 60.35011
29.45179 -8.654119 17.92464 49.99797
59.55883 4.149647 19.99140 66.80787
58.84823 -0.2417098 21.09238 65.55506
52.50514 -7.060569 20.63996 55.19292
50.28451 19.2696 17.77059 61.71615
59.03877 -4.343245 27.10155 59.8859
45.43183 -6.863089 12.28885 49.91697
49.32035 8.420879 20.89393 50.84106
50.13844 0.4999864 20.33221 41.75904
39.34384 3.193330 14.65490 43.1003
44.83127 3.980889 16.44049 48.04685
50.40867 -1.390823 20.30198 56.67151
36.19917 -12.26468 9.151291 46.72176
40.22976 8.719763 11.43870 41.32944
57.89613 3.30016 22.88141 50.66172
26.59688 1.680157 8.066935 23.10448
48.97289 10.49208 16.55740 49.23506
65.35547 -0.7951255 22.27252 69.34921
44.57912 1.899662 18.04871 51.01605
44.39202 -1.043773 16.9903 28.02483
54.586 0.3649564 24.53672 39.79551
54.43304 -0.1019753 23.70553 59.70726
49.14259 4.980435 14.78589 51.89241
45.76488 6.454239 17.76859 39.93045
52.95344 5.194063 17.22734 56.85135
36.44763 8.320695 13.95009 23.96988
56.23278 -6.391832 20.77251 50.95284
47.57699 4.790073 18.87454 54.58012
49.15628 8.409388 18.92746 39.71056
53.4933 7.252423 17.76128 58.50884
66.68983 9.754585 25.24448 76.18666
37.83893 -4.349798 16.15573 37.35371
71.97887 -6.061746 31.60432 52.605
43.87705 13.56080 22.47352 53.30198
26.76345 7.272119 10.96189 23.56650
33.43717 -6.120423 15.05788 27.5597
54.09092 10.02747 17.24765 56.39159
46.15204 0.5290926 21.17365 61.65141
55.24307 2.852504 27.11708 62.0633
50.84297 -4.294742 16.32535 54.61742
65.624 2.859510 26.82204 59.09405
65.85077 7.121138 26.74387 36.86199
43.98383 -9.29969 15.05746 34.56236
61.07763 -1.944802 26.23592 53.78754
36.19499 1.084941 11.50096 31.10465
47.69873 12.69054 22.54052 32.58539
52.44442 4.190038 21.12674 58.59621
59.88517 -1.250639 19.46250 65.7641
32.94451 -1.136390 19.51041 40.92487
57.79564 -8.407072 23.83641 61.67233
57.64336 5.609772 23.53299 27.25373
24.14377 20.84283 6.286495 17.65672
54.91191 -1.058408 17.18048 57.13659
48.01681 -1.416261 13.72837 45.51629
51.46387 2.482261 21.45468 75.73574
55.03506 13.23371 18.78182 62.12609
53.01681 7.32088 25.06812 53.46893
50.12581 8.700351 20.06746 44.97009
42.79887 13.81924 13.83654 48.06722
51.99636 -6.128464 23.92265 56.47989
29.73661 2.997832 8.150694 31.56814
52.51846 1.368375 23.94166 47.62229
40.32329 16.35087 12.75677 33.92927
46.23133 11.40794 18.62626 34.67508
49.23698 5.287289 18.02653 63.6912
63.657 12.13615 22.2143 70.08704
53.59651 -5.657185 18.11096 57.65388
42.57572 13.22200 21.33085 63.08201
51.24406 6.51087 24.10684 71.7702
32.31901 8.681594 12.41364 42.05136
54.00084 -3.435386 21.15679 46.28186
45.06625 -4.550324 13.73160 43.57226
48.55997 5.725301 11.91107 34.12021
48.26821 -11.80383 18.32051 35.82252
53.7636 10.19452 20.57763 40.48895
48.74251 -0.8477363 17.73027 40.85583
68.21562 3.792301 25.74067 66.89121
34.12522 7.830006 18.96062 45.32195
37.90025 8.684584 14.13703 36.98629
53.00591 2.106507 15.03131 62.31615
57.64051 0.2994314 26.83368 84.82863
54.10367 -6.483067 21.15202 64.51115
48.3755 8.44393 21.18105 37.49787
51.85371 1.043642 19.25040 60.9278
37.88257 -1.559072 13.68255 40.58769
58.34782 -4.399514 21.42981 69.36237
55.55458 0.3853118 23.01027 68.99905
45.01997 3.014769 19.11328 35.70126
47.4252 5.78935 24.91721 69.38955
38.33967 3.215241 10.13793 40.45346
53.89615 4.863428 21.74116 45.83579
34.86415 10.10497 12.8474 31.22481
60.92325 0.5722536 22.26861 59.98841
59.88477 4.908409 19.65403 54.43282
54.1092 15.86793 18.93589 49.63611
43.40765 -0.3759017 15.42585 28.58066
55.22262 7.930147 17.92699 47.73782
53.4557 7.246561 27.43898 43.97807
64.83686 -5.252063 27.14685 63.21693
56.11187 -3.194705 21.08142 38.69288
38.85002 1.824345 20.98474 38.87096
51.66953 7.324016 23.63456 42.91175
66.42632 8.658774 29.22023 72.33342
43.04705 19.48860 25.86647 55.39432
56.01276 14.57870 22.03388 53.55534
47.12277 3.508823 21.90701 53.55116
56.74827 3.794125 20.60873 68.87886
45.77865 12.46980 21.70479 36.5648
46.62047 4.502308 11.93615 54.10634
41.05144 5.667886 12.01179 41.78962
53.2453 -3.868202 16.93013 50.27663
53.7888 11.31978 17.11798 56.03127
30.87982 3.827967 11.30715 39.13347
61.44193 3.093164 21.16498 55.60321
40.83449 -9.376403 15.02020 32.62954
40.53524 -2.643674 15.39708 43.82684
46.46074 4.325226 18.57003 50.86629
59.10701 -1.705923 24.65460 41.54565
26.08755 -10.89046 10.21728 40.6581
36.78553 -0.3661094 11.48415 39.59677
59.38142 5.733641 24.16569 52.08247
73.46626 5.271732 33.34848 63.0421
70.30836 6.975592 22.51193 80.17147
26.81629 4.455363 6.256196 17.25833
49.73607 6.214709 17.68276 65.0581
44.8094 11.15783 17.58416 63.2894
64.4734 1.031661 27.53044 51.93996
33.88039 3.063151 14.07288 36.74695
68.85505 8.0014 33.05484 52.5104
55.39552 5.217852 18.97744 55.57307
56.68197 8.77524 17.21217 46.21657
55.15803 2.375904 24.63461 39.15838
49.37972 11.82569 19.43332 56.72875
54.05996 7.678194 18.93528 74.90348
59.62563 -0.5811889 24.24420 47.90656
39.09769 9.140207 18.05274 46.39661
50.11329 -2.455696 23.84894 49.50384
44.33938 8.814076 11.29521 60.16206
36.01459 -0.1962692 16.51724 31.37381
53.73882 17.98346 21.23376 70.16935
48.24547 1.707999 19.53287 68.1842
48.48232 -21.53443 19.54494 63.04124
55.12126 15.67938 21.94086 60.66847
58.10177 0.9693898 22.39573 62.83706
62.81068 -3.831048 27.98906 37.49122
54.77863 -8.951692 22.00743 59.86657
56.29413 16.97813 25.68932 54.08021
44.35275 -5.856065 15.18683 39.55388
42.26324 5.311419 11.32615 39.53598
41.77416 3.349529 13.46924 38.3943
59.07119 -6.546123 25.18052 42.87593
52.62234 8.706852 23.26501 53.91848
60.74666 10.30069 19.97020 59.62361
51.7901 6.787966 12.43232 50.54529
49.48841 -9.041835 19.91347 42.82379
33.21496 6.379343 18.47109 22.22286
42.83608 -0.3887477 19.26198 51.98982
58.90566 6.479899 22.34672 64.79377
38.86708 7.121925 9.412589 25.39843
39.06138 -10.18705 13.68964 52.73329
57.58497 5.286579 21.49572 61.26554
48.86454 5.683303 23.08288 53.09138
49.71984 3.228688 10.94475 69.37639
54.0442 8.549081 17.55736 54.79792
43.68507 15.31966 9.38958 62.06478
35.43669 1.366480 12.70296 31.77277
49.96935 6.204912 20.71118 53.70317
48.31328 -2.892094 19.08218 60.34956
55.42967 -6.080345 30.24604 52.66882
39.36169 -0.891561 14.22630 23.4967
44.73431 2.609196 15.91839 40.76198
42.34211 -1.385335 15.74328 47.39508
38.64424 -7.679498 9.553005 41.8961
49.37423 -9.645575 25.39812 48.2583
39.06988 6.532491 19.64858 47.28758
34.66840 -0.09675347 16.09624 44.98962
49.28633 15.36355 14.90732 39.89759
54.09937 -3.947975 22.33034 60.89902
38.98454 -6.131524 15.88354 48.90043
40.85694 7.367927 16.40521 25.05449
64.314 16.26035 24.28294 64.15392
48.57174 13.21193 21.93154 45.1065
47.8432 9.845876 15.83042 45.97286
56.27768 -2.011703 23.30932 61.55743
61.77799 -4.914055 25.48358 71.37419
48.40856 -1.979883 23.69773 45.7072
41.78199 6.615775 14.95757 27.73166
64.90359 -1.004709 24.43385 65.57712
63.58916 0.1358547 20.31097 49.51372
43.64577 5.670877 17.17293 45.57197
59.09745 4.490011 29.98890 47.95378
54.16604 10.40868 20.26598 56.02143
61.2897 5.421311 25.50053 46.19753
50.16458 2.223912 15.60141 55.70976
38.39669 0.9106178 17.42345 38.76754
54.68533 -4.205691 16.40946 65.26859
40.68546 -8.52779 14.70049 51.7189
44.77031 1.154359 21.65772 42.99853
55.78028 2.272698 17.14639 47.61636
52.97571 1.170537 23.05819 62.51508
44.53947 0.3084864 15.57685 42.11613
45.80788 -6.598917 16.28701 43.42877
47.33106 3.977682 15.7602 32.18884
81.2048 8.183425 29.92683 56.03272
50.76404 -3.066898 16.86412 53.28929
51.53195 15.04570 22.14058 52.74327
56.91241 -6.881064 21.32369 68.32285
44.02488 12.823 18.34239 21.94710
57.45945 -8.374578 23.02481 44.50366
48.42962 2.161985 21.95571 41.16837
53.17171 -11.51410 18.88203 51.80852
45.35711 10.79081 20.38162 29.50265
41.59632 -2.366603 17.23554 44.3918
45.37332 3.255886 16.42494 26.01951
66.28513 8.379297 27.27562 39.57251
55.04822 7.380156 23.33530 64.04498
43.72296 2.194734 17.80366 38.10659
48.65726 10.04495 17.72594 41.35633
58.93529 4.307534 21.44280 61.99892
57.77978 -1.462292 30.6931 67.02023
58.74133 -4.588407 16.35302 60.57455
25.32371 10.03417 10.42988 31.74891
57.62253 14.46439 18.87468 53.8298
70.71997 -3.034144 28.44621 58.7478
49.49285 -4.353228 18.58245 53.7977
46.48931 1.519150 16.62045 49.7155
53.54585 -3.508772 19.97289 34.74634
37.39493 -4.058478 12.00676 49.65412
68.68831 8.403163 20.04808 69.88225
53.95201 -0.8036485 25.7333 41.62982
43.82498 0.8353296 17.17056 55.96363
53.26219 17.63917 18.13693 58.99458
53.31943 12.46425 27.6435 59.8498
53.62114 -9.850218 16.39236 56.25944
44.22067 3.050159 13.34755 42.88112
54.0705 -14.17018 21.92201 68.51391
61.8496 10.48572 27.12702 48.5597
47.71082 10.56547 18.45355 48.03168
42.93826 -4.756591 15.26377 36.96124
57.00923 14.07864 24.9488 46.32931
47.52025 9.11744 18.42165 33.73318
60.7141 8.871266 28.66087 54.86872
24.47839 -2.786914 11.78060 31.82981
50.74243 0.6677095 15.15231 60.84145
41.77558 8.677417 11.28836 52.63567
63.34777 12.40588 23.62624 64.5393
62.01297 10.02279 29.36905 52.33118
61.11495 15.09229 28.13999 46.82123
41.87226 10.85900 12.55758 43.09227
45.57343 -4.951183 17.26845 37.31132
42.54236 6.836343 18.29875 59.27049
83.40856 11.03458 29.60942 85.00182
44.5255 -4.328198 19.70225 46.80362
59.01798 1.753738 22.32026 64.33994
56.09829 -1.235541 24.03586 52.40176
58.01278 13.71630 22.97000 75.98166
57.36326 4.591196 25.90256 53.21546
55.42411 -1.143736 20.78772 57.60663
49.93332 4.294792 20.43059 43.95717
38.12042 -3.429895 12.5381 37.35702
47.11148 -0.1824412 19.04586 44.34437
44.38587 0.5918665 23.66685 31.74505
39.18866 -6.614968 16.7935 28.35371
50.19105 12.86746 17.98825 40.86316
61.51687 2.97914 31.17259 43.78369
35.04929 3.929602 11.20934 41.74784
66.88222 -4.962832 26.49906 70.18137
58.62929 1.040376 21.04542 58.62515
44.36823 -3.847925 15.47678 28.41353
50.9544 3.411809 25.41637 53.79733
31.02899 17.29386 8.095708 18.90332
49.70619 2.645354 21.17133 52.21461
60.21066 -5.429708 20.25886 49.67136
37.78918 11.44529 11.62640 42.38456
59.59195 7.32593 27.72493 67.4737
59.48052 17.6721 21.47732 52.21488
55.59286 -7.043785 20.47393 57.11539
58.11713 4.091656 15.55736 43.8206
59.74675 -2.773232 16.31989 55.5046
63.23438 -1.702270 29.81902 52.65655
41.07868 7.789948 21.30436 33.57119
41.06975 8.617235 20.94633 46.8047
47.26647 -2.736527 17.51564 47.98941
33.76289 -0.4600734 14.06749 38.29962
58.49543 8.540616 20.92075 58.55932
46.17608 -3.332358 19.84594 39.76436
55.92106 4.650755 16.06241 55.87971
36.68021 8.695798 12.82582 34.88202
54.44295 9.37359 22.63770 57.18928
39.85584 -0.3960061 7.079331 41.09933
46.81241 7.222772 16.75245 58.42668
61.71044 7.235608 18.19497 66.81891
69.99467 5.67094 27.374 50.99303
60.33253 -1.703566 27.72115 58.50857
47.24514 0.8480002 16.83138 46.12023
49.37226 21.64616 15.56421 42.69047
39.31795 -3.939138 12.88887 37.94769
56.34995 -2.604590 24.78852 59.01998
62.99316 15.11994 25.12581 58.86634
27.67895 -2.646529 11.61963 22.42947
41.6029 6.918429 10.38468 40.02674
55.05527 2.373012 21.48015 50.85709
35.63742 -2.450121 14.98572 22.75035
55.37802 7.22383 25.15293 74.67586
55.74853 -1.569095 17.72058 52.86253
33.14594 13.43592 13.49071 19.86052
53.97746 -3.969389 19.45272 69.90991
42.6752 20.88976 12.99029 42.19767
40.67434 -5.279338 17.22070 46.01282
44.85665 14.05799 13.13810 50.20145
47.39532 4.850963 15.52587 50.23565
44.75813 0.4860118 17.05484 49.42128
43.63428 9.317992 11.17524 54.96257
51.63514 7.12898 16.90574 76.3345
44.892 -2.373124 16.69128 44.14674
49.45191 -8.287656 19.30231 54.90011
61.62966 7.274818 18.89314 69.5845
52.80673 15.42787 23.57102 32.1874
52.03011 -0.759945 19.02679 45.1055
36.28238 -7.329426 15.19359 50.96374
53.79494 0.6901786 19.26644 42.53443
37.36745 0.6343899 13.89918 39.68382
50.85536 7.803232 21.81077 38.89775
44.6144 -0.05662501 15.21081 24.02865
45.15917 10.76237 11.81272 44.63008
52.2677 10.26474 14.66819 59.63404
45.4052 9.164518 15.8748 27.49622
50.09451 9.759 14.63492 25.90424
41.21958 7.620678 12.21550 54.3206
49.07768 7.14001 19.39317 52.59705
51.22407 0.4873102 19.98736 64.3181
69.89465 4.284478 23.13207 71.468
65.9434 0.4386691 23.91611 81.6142
53.29608 13.14064 14.32771 41.25041
52.67964 8.398883 19.33138 49.70049
35.71905 -4.690741 12.90544 37.13983
65.31073 5.980597 21.62468 57.607
49.77643 -4.705812 23.74905 57.24223
36.60589 3.162954 15.39453 28.63667
49.99168 1.376404 16.47252 32.57689
55.0787 4.232147 19.87869 54.42524
41.81152 2.255201 18.28535 24.96889
54.28878 -0.9765505 21.62945 45.58172
50.79173 9.993296 21.16342 57.21549
46.12364 -1.532176 15.43658 55.36385
55.133 5.518339 23.63813 68.88534
47.80002 -1.191592 20.43180 52.33085
42.45264 15.84223 13.95406 44.85478
58.27121 13.67436 21.12923 68.07198
49.11102 7.95799 20.81707 52.3629
51.16691 3.037852 14.73927 47.58764
41.98612 -1.079454 13.57764 38.76133
64.41549 4.157859 23.72416 45.40483
29.9604 4.645774 14.68083 27.3442
51.02253 -2.508397 16.36751 49.74401
38.8691 -0.2028261 12.86854 35.64026
56.9565 -3.203786 21.32083 68.25416
37.53949 -5.298522 14.05683 33.79837
48.51417 -8.584027 17.28389 56.87966
57.27387 4.545264 26.01884 68.3199
40.9367 8.388785 15.04107 43.22461
38.86728 6.795023 14.89140 38.71906
50.62005 -0.928655 22.43378 50.98833
55.73224 3.8226 22.60302 47.98252
49.337 12.36312 21.44033 49.09599
40.94014 6.019235 22.63988 27.37827
54.82275 17.14067 21.22033 60.40025
30.26392 0.06316697 11.52665 43.21693
54.54599 8.85506 19.5567 66.01324
47.69885 -3.82504 23.47038 46.65729
53.65547 19.41216 27.81673 47.30365
56.15439 3.681035 24.91513 68.69313
75.22148 11.73966 29.04425 62.73034
40.38708 1.203400 20.48614 51.68565
42.54723 3.818939 22.78269 33.76868
67.13459 -2.134234 21.72449 66.22025
45.68465 9.342554 16.79275 42.22965
51.68726 -20.31493 20.20590 52.49826
60.66892 -9.236747 21.58120 65.21285
56.13715 0.1408979 26.09226 60.7464
54.49441 3.272842 22.31331 63.38593
47.1787 -0.614443 18.98358 45.23967
61.01369 10.53956 22.03796 61.15634
33.0434 1.710027 8.318456 21.18222
48.01317 -4.890684 12.36127 52.82928
54.42551 6.749779 19.09129 57.064
55.25833 9.9348 27.09836 43.54478
67.28864 2.290348 27.45477 63.45227
59.81277 0.1539041 22.50517 62.71646
43.78275 -3.20156 17.03778 27.75368
65.13317 -0.5028082 29.53686 62.15449
54.68977 5.426791 28.44692 55.1669
54.33469 10.86265 15.80755 52.85767
43.86489 0.9777324 18.27711 47.75802
59.88471 12.80257 24.9312 72.6737
55.04461 9.92596 15.34458 65.44657
47.87519 0.7445574 12.32965 72.05105
42.64161 9.631529 17.66188 31.70682
44.03236 13.2341 20.78058 44.44485
56.80818 8.498717 27.27883 51.76617
47.31621 7.516062 19.04743 34.08694
34.0247 9.601698 14.06770 29.22559
49.22712 2.890787 19.16268 46.1354
50.72603 -14.20376 20.99923 64.7304
48.80112 -1.335916 15.04928 51.3007
44.88838 -0.94571 15.48166 52.78398
30.12507 -9.205157 14.35123 33.56024
51.86355 -1.072048 21.99265 56.83661
63.15042 11.25744 22.84697 63.15531
57.10962 6.578032 25.95038 62.91333
50.7127 7.944212 19.82808 66.4457
38.84542 -1.099172 9.806064 49.62192
55.85912 -7.285945 26.31835 40.81715
40.69722 8.076507 18.90315 32.63169
51.42381 10.20819 17.67442 64.21663
63.50688 1.948752 29.31103 48.55231
46.51957 -0.9314324 22.28629 46.55982
54.65897 0.2908350 19.9915 62.47015
54.32166 13.64623 21.73247 62.89811
38.73943 -1.190599 13.5414 30.04813
53.91912 7.169232 17.14703 69.34605
47.71163 9.813518 17.85789 53.84949
62.72772 15.30444 24.59063 55.9661
40.30492 -1.178631 14.86803 44.09741
46.39432 -1.057279 19.11413 40.49707
56.74309 18.87834 24.45470 56.46179
51.79315 -2.946906 20.83218 60.79951
42.10143 -3.789814 13.83079 49.51949
46.17645 6.531711 14.71591 59.79538
52.69893 -3.340308 20.16043 53.86045
48.68806 1.754878 17.12362 33.3709
41.92036 -0.4378864 18.96018 57.20926
47.33391 9.903861 19.94968 46.24269
75.41595 12.77358 26.91113 83.9159
56.92331 12.35901 23.51364 34.29779
41.36600 4.092266 13.67126 55.34081
58.95643 4.905904 20.30149 57.05644
46.48007 8.40095 17.98870 43.50009
69.81733 7.626208 23.85312 67.29803
59.51663 2.098946 22.32835 41.23403
29.88096 -11.11723 12.16694 37.8223
53.35496 -1.763543 22.04936 62.47231
47.67394 8.733833 15.49805 37.2636
53.12558 4.807173 17.16524 66.888
57.35365 1.20005 22.93623 67.79992
60.18102 0.3326909 22.30841 59.54785
66.04277 4.975833 24.59072 73.9767
52.08193 5.864035 17.42318 46.63601
58.3072 -4.283239 20.45064 51.57332
51.31144 -12.23417 20.29298 48.63341
48.4008 3.61443 20.41035 38.83863
57.0495 8.693174 15.93369 58.34179
41.78283 8.504102 14.63974 44.6281
26.37981 -9.5862 5.902373 50.98046
38.2025 -1.436279 10.75539 33.93151
39.02257 -0.0795461 13.68205 30.90802
63.78602 -6.229248 23.37354 50.33398
41.6225 6.347579 18.83049 60.01397
50.98083 17.61319 20.7728 49.17911
61.12413 9.138548 17.32693 61.62614
25.53722 4.934389 3.60641 35.38249
65.87488 15.60881 22.14230 75.30191
47.47258 11.2685 21.535 22.7381
52.13593 -1.391562 22.05756 59.25493
42.85993 -7.028862 17.2528 61.70683
52.76418 5.881355 19.16676 46.83142
44.27458 4.474306 15.36534 48.64894
46.75474 6.111763 22.39003 49.31481
40.03725 7.703231 12.63798 31.04878
35.28821 -4.876697 17.06257 65.49412
53.20997 6.899141 11.95647 50.38232
39.63746 -11.89098 5.782316 42.34723
48.92189 12.73269 25.75742 54.72582
41.85545 1.452485 15.57795 40.88776
54.78765 7.539781 24.64020 48.02022
57.59353 7.216395 26.89426 66.16208
47.219 0.98581 14.63727 47.92469
51.20543 0.09259385 15.40908 43.74882
45.57773 10.27622 12.05995 52.76741
50.94095 8.593899 21.52522 33.35471
49.41009 5.033688 19.87464 68.00966
32.99206 -8.25807 10.89571 31.7013
43.00241 -4.55001 16.82034 48.2274
55.85164 16.50756 21.21968 65.50456
51.94505 4.311148 18.58981 49.51577
22.29264 10.19321 8.807824 17.80639
41.29039 7.603837 21.31213 24.72423
40.98491 7.232681 18.15716 42.69361
39.00043 -2.448485 14.81503 47.57224
46.90163 -2.104624 20.51735 30.61411
34.4148 7.107138 17.40879 33.62787
36.39398 -3.095294 15.46363 36.31153
53.75927 -0.8825684 15.44426 57.34617
40.46641 -7.851308 15.35879 33.87549
48.54501 2.895400 19.40896 56.74133
53.51811 -3.112738 17.23571 49.37379
49.89715 0.3876615 23.03248 54.27001
48.15766 -8.721644 22.32676 44.30458
46.16639 1.930796 12.63313 39.61946
50.94754 4.243396 18.31385 48.84745
36.86526 4.931829 14.48334 26.91154
50.64524 4.113626 25.08682 49.00419
28.74042 -2.871925 11.04577 10.79332
55.9087 -3.210486 23.07479 37.80772
37.73887 13.87986 15.71202 41.62297
24.32980 18.11231 12.86661 19.73414
55.39625 -4.763357 19.99735 59.79372
57.16707 3.686038 26.60973 45.47296
59.51321 17.39736 26.13501 70.75698
56.36915 5.823617 27.60674 45.71415
43.91938 5.613671 22.95167 50.78551
42.85471 7.34328 13.44399 42.59938
47.62828 -0.2286385 16.95952 38.35376
68.75052 14.10395 22.38109 82.23013
56.62596 -0.676623 22.39061 52.78448
49.23282 7.187152 15.46581 49.49732
35.56163 4.864949 11.66103 47.30657
47.64718 4.224139 16.92794 53.13547
57.24023 9.693093 22.14284 64.47852
40.09588 10.44241 14.84347 45.89301
48.19389 8.240209 24.54570 55.06081
35.06633 -12.30828 13.05675 55.79669
35.51001 4.968912 17.18222 45.78132
49.21156 0.6832834 15.44924 48.5414
66.56935 -6.973784 24.53541 70.08622
38.04989 -5.532575 15.47116 37.87041
55.67546 8.600036 26.54975 48.49944
39.95599 7.040959 14.66411 57.70497
43.48450 9.456422 22.64099 39.60096
33.82183 7.50383 9.532232 32.70466
42.85369 5.69953 17.53943 43.57636
54.42622 8.795007 19.97518 44.71488
61.56923 15.60200 19.54735 68.26665
44.93665 5.365647 15.3018 50.87753
50.61438 14.54369 23.40480 51.89041
44.48844 3.742543 13.78786 63.1765
58.4448 -11.29651 22.64451 45.71587
38.52331 -15.74596 13.38931 39.8719
51.11058 -6.131575 19.71275 38.27245
59.0543 4.573761 22.39492 33.22500
39.45927 -8.397676 15.39175 34.49714
49.0542 -1.899996 14.94644 34.30608
52.96821 0.8465435 24.00737 59.76748
37.31856 8.498543 13.68324 34.73792
42.86177 8.966008 21.06236 25.42074
46.94618 11.26505 18.59283 49.20128
30.25900 -9.73547 12.67027 17.03394
38.62156 -4.717086 15.84583 41.32916
61.35311 -4.765051 26.90832 51.92822
49.59635 3.047772 23.48041 51.40419
52.02974 2.293412 20.37837 67.17684
50.19866 2.375115 19.91393 57.71842
46.03202 -4.013739 17.67725 61.51833
56.86154 14.15691 22.48793 43.30352
56.41994 -7.764509 19.64890 53.54802
54.06958 -3.601526 25.26686 59.35605
56.83714 -9.573361 23.74937 57.70736
53.64752 2.096009 20.35184 57.93559
37.27191 15.30242 20.21106 33.38579
59.88905 11.07137 17.94779 71.80114
55.42092 -9.012312 18.61974 54.39239
39.05826 -2.005509 14.99925 54.7657
39.92606 6.116922 17.55380 15.37107
35.33501 6.035102 18.93004 33.78016
43.38075 -0.8498992 11.55096 31.35522
51.15989 3.137947 17.99312 60.41017
55.67945 -3.199316 26.5845 68.93158
44.06857 4.256001 16.96168 55.44931
72.77981 21.3873 29.13402 63.98076
62.06794 -2.856957 23.64693 58.15178
65.27273 9.363268 27.89871 79.48394
55.48338 1.583964 28.82452 60.35442
49.26494 3.065723 23.86022 35.78034
37.16467 -3.88511 9.301194 32.40567
54.7703 15.03396 21.67751 37.43552
58.50714 13.30801 19.32692 67.2622
40.02551 10.48335 9.715556 21.38574
46.17735 11.74318 16.51314 45.3918
50.65409 -0.3688263 15.42613 43.45416
36.85831 7.61677 12.41689 40.44595
65.48127 15.54783 30.48558 60.1603
43.28649 -4.36009 16.61700 40.53279
47.33055 -9.820944 22.52366 45.69401
60.24609 15.18745 22.2615 80.29985
47.6109 6.291483 18.42406 46.20005
32.52579 2.163241 9.673215 35.27061
48.10742 5.023338 28.71845 44.01503
42.53343 10.37183 22.23493 51.5103
51.65264 0.7802082 15.98401 53.49945
73.92167 7.71875 29.47882 65.73937
50.08472 7.385094 14.87745 50.5794
31.05070 5.02159 6.320998 47.83895
50.46649 3.390473 17.70616 48.57303
65.88473 15.99298 27.95033 48.00119
71.943 20.16625 26.64881 93.15871
50.15365 -5.470214 21.19497 51.69552
53.1029 4.348359 20.93293 25.96702
57.25345 5.215286 19.01661 43.09172
47.21926 -18.31202 19.56487 44.35467
51.86111 -2.302820 23.37832 34.78882
40.86635 3.55813 11.12562 59.46731
41.81289 1.341042 8.254875 53.02104
66.34576 4.365926 29.19620 53.34384
53.60574 -2.281706 25.03334 53.4722
53.36986 5.65597 15.21477 60.28249
44.3015 10.90607 16.77324 63.7355
42.04445 1.970218 22.04973 44.00357
52.80729 0.551982 15.49625 54.89782
59.94559 9.41296 19.90171 45.41671
55.74521 16.31075 26.3625 60.29774
37.11364 1.525146 14.04184 34.33557
32.83712 2.774706 6.601578 41.34569
43.33734 -3.330938 19.47048 38.48738
47.11865 8.52896 23.66342 62.64756
59.86298 0.1594297 26.43365 72.1182
47.38146 11.47643 17.11502 56.6119
50.51492 -6.982039 21.46275 57.37549
41.04778 2.430464 11.45823 31.58703
47.56391 8.033224 18.18477 52.83157
38.83323 5.135054 18.99735 33.14404
33.55940 4.576157 11.44899 28.96809
65.17481 9.532619 27.40959 53.81124
60.97809 10.72855 23.05310 49.10621
65.71203 4.680273 22.60413 65.91371
44.6671 12.72787 15.06389 37.87554
46.70985 4.207429 15.65107 49.61584
44.0126 -0.8969562 24.53429 37.99542
56.09004 6.60307 19.38070 62.44445
52.74047 2.275871 22.03446 48.78008
43.1794 -1.274556 18.38257 52.86725
52.52871 4.785603 19.78508 55.5101
50.91668 8.13449 20.70597 44.02226
53.04344 22.67535 21.92497 42.16448
54.31046 -5.32972 21.98447 64.63553
51.76945 4.045078 20.51839 52.57333
57.46368 4.722063 26.79969 67.60241
48.26568 0.9973341 19.01231 43.01635
39.52084 3.624995 20.03794 19.2539
56.47525 15.15424 18.85206 54.79408
64.2143 7.907288 25.63994 69.30486
46.54991 -5.522777 21.88838 50.19121
48.57942 9.917584 5.490708 58.09242
68.94336 15.78718 32.37786 71.80714
65.097 5.148921 30.7213 84.25914
55.61453 2.618225 14.70006 58.49522
43.70043 2.583815 21.14595 39.04892
62.93275 17.63986 25.88950 75.09722
40.43976 -0.585903 13.92127 36.75159
55.45437 -0.3992024 24.08535 36.08245
43.35162 9.897227 14.66669 42.22049
47.67736 -0.7712326 15.27691 45.10154
33.92403 8.630114 16.55974 31.888
54.21994 2.466792 20.30346 55.51824
44.68194 8.934485 15.75736 23.64034
51.26517 1.499947 22.02061 39.83077
37.44306 -5.439318 8.345963 38.25729
44.72494 5.445372 22.24458 67.62202
50.8323 8.911907 20.07359 57.16967
55.58463 1.375984 24.00998 58.9103
51.32851 12.48013 14.64975 36.59635
50.01345 3.827587 20.03792 49.61945
39.40343 14.19370 15.54872 45.10012
44.0108 -6.723932 14.32292 61.41351
44.85873 10.95146 16.78062 51.02501
45.62154 13.41332 20.0166 47.36642
46.49261 3.184934 21.42679 63.42772
21.25775 10.71705 10.84546 33.68133
72.93982 8.190522 27.71297 75.4493
65.63378 5.692074 27.75515 55.73338
47.33759 3.022367 13.15211 44.59232
38.77872 6.00241 19.44568 48.34907
63.15298 -3.808 30.82053 80.11868
43.96191 -3.850284 19.4795 43.23053
44.00066 9.85134 20.63579 55.01284
50.88649 14.63406 23.84008 44.95536
50.22891 -4.04665 20.18481 32.5259
43.60246 -1.379428 16.41767 40.18411
48.53157 5.588242 18.69031 61.24137
41.77733 1.084981 19.79319 43.38450
58.40972 8.500872 24.77230 49.15017
57.57225 -3.797773 22.58461 69.06964
40.20747 -4.917563 13.88889 47.72803
51.50455 -9.029045 18.24453 46.40172
51.22361 0.7188372 19.85480 44.53326
42.90079 10.53299 18.18231 35.82266
54.97448 6.566912 21.16093 38.28001
44.30551 18.52872 16.14252 29.35241
42.90053 11.1345 20.60380 47.55508
40.59100 3.836315 19.17098 40.01481
64.46993 3.410776 19.69802 52.40692
46.05315 3.149181 15.95172 37.64019
38.9236 13.07707 16.62506 38.5357
50.07179 -2.132767 10.71973 53.55615
62.8026 3.477998 22.30143 61.35994
29.33549 -2.266676 8.703178 28.27489
52.82056 -5.803662 20.60624 63.00311
39.763 8.74093 11.52920 48.11935
49.20929 10.07257 18.14720 47.40108
59.08706 12.83303 28.72768 63.30883
48.69612 5.591273 22.10258 37.92279
45.22457 7.19037 18.16479 29.70533
58.29379 7.562675 21.38561 56.18841
60.23416 0.4790846 25.61319 54.84405
40.80870 14.72674 16.83701 52.38186
44.07144 5.012861 21.17531 41.13701
41.01263 -6.398877 21.61346 39.50615
57.26699 12.37304 21.68357 66.50561
45.10022 3.784656 17.50116 46.27478
47.70641 -1.740854 21.81915 49.36876
49.49623 7.38965 13.35020 53.29553
67.39344 0.3582493 32.50292 83.45366
43.91011 6.955515 20.1076 42.81388
51.06955 11.46009 27.47350 61.35225
70.06126 11.46235 27.73504 85.74356
62.28284 7.52419 30.89928 64.75573
46.35349 12.13339 13.76977 37.60365
61.63559 4.779976 28.43432 80.95307
49.11642 7.717767 23.40764 62.98167
29.38953 15.28556 8.171793 24.36978
59.50431 0.1827939 24.80716 54.25433
47.56914 -5.181217 17.58222 43.88487
49.92213 7.02148 15.47130 61.99712
49.03364 4.933716 12.89841 39.03318
45.31894 8.319285 18.97337 47.01734
46.90397 10.42634 17.42944 48.40244
38.53745 -3.053576 12.14128 38.62005
53.53033 -13.54362 25.11791 49.2909
43.88261 -3.446498 16.81401 51.68662
56.41206 10.91076 23.93817 59.53211
50.15494 9.85292 18.61940 53.80352
53.84678 -1.913277 20.18228 43.39612
54.07992 8.641245 20.06299 44.45265
47.44463 2.684312 18.93508 48.09456
47.47627 10.80877 17.72176 58.99301
39.82833 8.556109 14.18629 55.5331
33.24262 -4.679348 13.79564 20.6876
53.47825 2.301145 19.61014 38.3461
43.7571 10.22414 15.54398 54.45097
39.99939 -2.742248 22.66244 37.48000
62.0932 2.095559 24.75525 69.32196
25.84546 2.086693 6.911238 34.54912
54.73739 -8.796734 15.48919 43.64274
46.72121 11.00066 15.6844 52.25666
47.94386 11.10268 18.31938 44.31964
43.60601 -6.907512 12.91284 45.1927
60.80403 -4.68798 22.13393 60.51706
42.9251 8.836826 16.09451 55.54047
41.19792 7.915696 20.71839 47.43137
39.85169 -4.470165 16.78824 50.02049
44.74882 6.133915 17.40697 45.12814
51.0531 1.828886 16.24041 64.59942
42.4202 3.107667 12.01330 52.86244
54.34774 -11.8652 23.02954 52.83223
35.47003 12.79844 17.7721 60.13769
49.39028 -13.65518 22.57524 43.38453
43.71669 -6.957278 17.17567 46.32105
58.35255 9.28393 23.99462 57.64893
76.23737 9.186005 30.78686 79.4134
46.33435 16.42037 14.99501 41.66136
51.12072 3.630508 14.11113 66.11836
34.31019 1.043082 14.72027 36.24138
51.97107 7.015492 14.57453 66.39966
62.43678 7.683385 26.39867 55.5293
57.27632 1.404660 22.24053 70.24073
49.5025 7.503132 22.38669 33.98815
45.59278 -1.347406 18.73675 44.78555
28.63468 0.8145263 14.24592 15.69798
35.79533 -12.44964 13.90434 28.42654
28.59595 7.324136 12.97275 32.70563
52.87595 -5.731458 23.34969 33.35862
43.8342 12.19212 15.74296 47.8095
50.24297 4.486603 15.69612 53.69657
50.53234 -11.37676 20.09458 37.95678
40.65528 -1.754301 21.72338 29.18788
50.68405 1.010595 18.3277 56.25837
39.82078 -2.874696 14.28346 28.24900
53.22794 2.402170 23.42319 51.68007

Read more!

Wednesday 1 August 2007

Statistical Modelling: The Bits Pt. 2

In my last post, I outlined the ideas of statistical inference: that we have a model, and we have a philosophical scheme that tells us what maths to do to estimate the parameters of the model. But this is not the same as doing this in practice. The equations are often complex, and difficult or impossible to solve (anyone fancy trying to integrate in 1000 dimensions?). Because of this, a lot of methods have been devised for estimating the parameters, and particularly estimating their variability (which, from now on, I will write about as estimating the standard error).

For maximum likelihood, the problem of getting the estimates is an optimisation problem: find the set of parameters that give the largest likelihood. For the simple linear regression, this is easy, because the equations were worked out a couple of hundred years ago, by some French mathematicians trying to work out where Paris was. For more complex models, such as generalised linear models, there are algorithms that have been shown to efficiently iterate towards the ML estimates. In more complex models, more general or complicated search algorithms might be needed (for example, there is the EM algorithm which is useful for dealing with missing data. So useful that the missing data is sometimes invented, in order to use it).

These algorithms only give point estimates. We also want to know the standard error, and there are methods for doing this. Again, it may be that the equations can be derived directly, or can be found with a simple numerical search. For less standard problems, more complex algorithms can be used. Two popular ones are called the bootstrap and the jackknife. They both work by re-sampling the data, and using that variation to estimate the standard error. What is less well known is that they can also be used to estimate the bias in the point estimates.

Even more complex problems have even more complex solutions. For Bayesians, the most common method for fitting their models is a technique called MCMC. This estimates the parameters (actually their full distribution) by simulation. The technique is not simple to explain, but it is used because it is very flexible and for a lot of problems works well in practice. However, this is not the same as a Bayesian analysis: MCMC can be used by frequentists as well, and other complicated methods can also be used to fit Bayesian models (e.g. sequential importance sampling).

So, now we have a model fitted. There is still one part of the process that is often done: model selection. This is not always done, and is perhaps done more often than it should. It also creates a lot of heated debates, which may surprise some of you. That's because model selection is often disguised as hypothesis testing.

A general description of the role of model selection is that it is a way of choosing between different mathematical models for the data, to find the best one. For frequentists there are two principal ways of doing this: null hypothesis statistical testing and using information criteria to compare models.

Null hypothesis statistical testing (NHST) is the method everyone is taught to test hypotheses. What one does is set up a null hypothesis, for example that there is no difference between two groups, and then test whether this hypothesis is supported by the data. If it is, then the null hypothesis is accepted (strictly, it's not rejected), otherwise it is rejected. In one common formulation, whenever a null hypothesis is set up, a more general alternative hypothesis is also proposed: rejecting the null hypothesis means accepting the alternative hypothesis.

Where is the model selection in this? We can take the Harry Potter example, and test the null hypothesis that the change in book length over time is zero. The relevant part of the mathematical model is

mi = a + b Xi

where mi is the expected number of pages, and Xi is the year of publication. The null hypothesis is that b=0. But this is equivalent to this model

mi = a

And mi = a + b Xi is the alternative model for the alternative hypothesis. Hence, hypothesis testing becomes a problem of comparing two models. The models are compared by calculating the amount of variation in the data explained by the two models, and seeing if the larger model (i.e. the one with the effect of time) explains a significantly greater amount of variation.

The alternative method for model comparison uses information criteria. These are measures of model adequacy: they consist of a measure of model fit (the deviance: the smaller the deviance, the closer the fitted model is to the data), and a penalisation term for complexity (the more parameters in the model, the higher this penalisation). The different models being compared are fitted to the data, and the relevant information criterion (e.g. AIC or BIC) is calculated for each one. The model with the lowest criterion is declared the best. If there are several models with similar criteria, they are sometimes all examined, and one of them chosen to be the best, for other reasons (e.g. because it makes sense scientifically).

Both of these methods have their counterparts in Bayesian analysis: hypothesis testing can be carried out using Bayes factors, for instance, and there is an information criterion for Bayesians (DIC. Alas, not the Bayesian information criterion!).

So, there we have it. OF course, this is not all there is in statistical inference: for example, I have not dealt with model checking (i.e. asking whether the model actually fits the data well), and there are many details I have left out. But I hope this give a scheme that separates out the different parts of fitting statistical models, and reduces some of the confusion. To summarise, here is a list of the parts of the process, and examples of the actual methods used in each part:


  1. The mathematical model


    • regression, ANOVA

    • ARIMA (time series)

    • Generalised linear models, and Generalised linear mixed models


  2. The inferential framework


    • Frequentist (Maximum Likelihood and REML)

    • Bayesian

    • least squares

    • non-parametric


  3. Model fitting (i.e. the computations)


    • bootstrap and jackknife

    • MCMC

    • importance sampling and sequential importance sampling


  4. Model selection


    • Significance tests

    • information methods (?IC: AIC, BIC, DIC, etc)



Read more!

Monday 30 July 2007

Statistical Modelling: The Bits Pt. 1

It's a sad fact of life that it's too short, so we don't have time to learn everything that we need to know. For most biologists (and probably most scientists, and indeed social scientists), one of the things they should learn to do science well is statistics. Of course, the poor dears are too busy learning about things like biochemistry or PCR, to have time to learn about the important stuff.

Because of this, there is a lot they don't know, and also a lot of confusion about the bits they do know. I therefore thought it was worth writing a few posts about some of these issues, to try and dispel at least some of the confusion.

Now we're below the fold, I should admit that I don't mind too much the lack of knowledge - there are large swathes of statistics that I don't know about either (I was too busy admiring the paper aeroplane building skills of my fellow biochemistry students). The confusion is more of a concern, and I think it is largely because biologists get little training in statistics as undergraduates, so don't have the grounding in the theory when they start doing research. And, yes, I am generalising - some biologists do become very good statisticians. But there is still plenty of confusion.

A lot of basic confusion comes from not understanding the different parts of using models to analyse data. For the moment I am ignoring things like model checking, and how to interpret the results. This is not because they are unimportant, but just because I want to write a blog entry, not a monograph. Anyway, There are, I think, 4 parts that can be separated out:

  1. the mathematical model
  2. the inferential framework
  3. model fitting (i.e. the computations)
  4. model selection

Point 1 is often passed over, and points 3 and 4 are often confused with point 2. I will deal with the first two points in this post, and treat the others later.

To help disentangle the different parts, it is useful to have an example. For this, we can use a simple linear regression, such as the one I used here. The mathematical model is self-explanatory. For the regression, it can be written as

mi = a + b Xi
Yi ~ N(mi, s2).

The first line is the equation for a straight line. Xi is the covariate (year, in the example). This is then related to an expected value, mi. The second line says that each data point, Yi, is normally distributed with mean mi and variance s2.

Note that the data are assumed to be drawn from some distribution. Hence, they are random. However, they depend on parameters, in this case the mean and variance. Here the mean is modelled further, using an equation, so it is deterministic. i.e. if we knew the parameters, we would know the true value of the mean.

More complex analyses also have a model underneath them: the usual ANOVA actually has a model that is almost the same. The mean can be a function of several parameters, and might not be the mean, but just a parameter that controls the likelihood. The variance might also be modelled. And the parameters that are modelled could themselves be functions of more parameters. And all these parameters might themselves be random. Yes, things can get complex.

The problem, then, is to estimate the parameters. There are several schemes for doing this. For the model above, the two commonly used schemes are the frequentist and the Bayesian approaches. Each of these gives us a way to estimate the parameters. They are, in essence, philosophical approaches to understanding what we mean by probability and randomness. In the frequentist scheme, what we observe is random, with fixed parameters, and hence any statistic we calculate (such as the slope of the straight line) is a random function of the data. In contrast, the Bayesian approach is to say that anything we are uncertain about is random. Hence, we treat the data as fixed (because we know it: we have observed it), and the parameters are random.

Both of these philosophical schemes are useful in that they then can be used to construct methods for estimating the parameters of the mathematical model. Both use the likelihood, i.e. the probability of the data given values of the parameters, albeit in different ways. For example, the frequentist scheme tells us that we should find the values of the parameters which give us the maximum value of the likelihood. Hence, this is often spoken about as a maximum likelihood approach (there is another philosophical scheme which leads to the same equations, but to me it always looks like asking Fisher to wave his fiducial wand to make everything mysteriously right. I assume I'm missing something important).

Estimating the parameters is not enough, though. We also want to know how reliable they are. It makes a difference to our estimate that for every year J.K. Rowling writes an extra 300 pages if the range of possible values is 290-310 or 0-600. In the former case, we have a pretty good idea how much time to book off work to read her latest book: in the latter, we have little idea about whether we need to have a headache or major heart surgery. We can summarise the reliability in a variety of ways: a popular one is the probability that the estimate is less than zero. Alternatives include giving a range of likely values (typically a range for which there is a 95% probability that the estimate is within it), or a statistic such as a standard deviation to summarise the variability in the estimate (in the frequentist scheme, this standard deviation is called the standard error).

In my next post, I will describe the model fitting and model selection, i.e. how these schemes are used in practice to make inferences. I suspect this will include a rant about p-values.

Before I stop to work myself up into a frenzy, one last point is worth mentioning. The approaches discussed above assume that the data are normally distributed. If we discard this assumption, we can use different approaches. One, called least squares, is based on optimising a property of the model and data (i.e. minimising the sum of the squares of the differences between the data and the expected values). It turns out that this is the same as maximising the likelihood, i.e. it is mathematically the same as the frequentist approach. I, therfore, use the frequentist interpretation of the equations, and it gets round a lot of fiddly problems. Another, called non-parametric statistics, does away with any distributional assumption, but also throws away part of the data: rather than use the data themselves, the analysis concentrates on using the ranks of the data (i.e whether each Yi is the smallest, second smallest etc.). In some cases, I think this can still have a likelihood interpretation, but I must admit that I haven't checked: it's just a hunch.
Read more!

Still here?


Wow! The virtual world is still here!

I'm co-teaching a course with Bill Shipley in (sort-of) central Finland. Yesterday we drove up here, and settled in. I tried to get a wireless connection (I had been told there was one) from my room. My laptop would search and find a weak connection, and then when I tried to connect, it would disappear. Repeatedly. Once it decided to mock me by giving me two networks that I couldn't connect to. It would tell me that they were secure networks, which is a hint that I might need a password. But, being rational, I tried to connect anyway.

This may not be relevant at all, but I've just been reading a book about happiness, which explained that we often do things that we expect will make us happier, but we don't get the thrill at the end. This leads to addictive behaviour: we just can't stop ourselves even if we are objectively aware of the futility.

OK, must stop: I have to check my email.
Read more!