From 48b978b6dd7774e754abc04529bfd9d041284587 Mon Sep 17 00:00:00 2001
From: Roger Bivand Spatial weights
library(terra)
-v_eire_ge1 <-vect(eire_ge1)
-SG <- rasterize(v_eire_ge1, rast(nrows=70, ncols=50, extent=ext(v_eire_ge1)), field="county")
+v_eire_ge1 <-vect(eire_ge1)
+SG <- rasterize(v_eire_ge1, rast(nrows=70, ncols=50, extent=ext(v_eire_ge1)), field="county")
library(rgrass)
grass_home <- "/home/rsb/topics/grass/g820/grass82"
initGRASS(grass_home, home=tempdir(), SG=SG, override=TRUE)
@@ -398,7 +398,7 @@
Spatial weightsinv_dlist <- lapply(dlist, function(x) 1/(x/1.609344))
combo_km <- lapply(1:length(inv_dlist), function(i) inv_dlist[[i]]*prop_borders[[i]])
combo_km_lw <- nb2listw(grass_borders$neighbours, glist=combo_km, style="B")
-summary(combo_km_lw)
+summary(red_lw_unstand)
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 25
@@ -434,7 +434,7 @@
Spatial weightsred_lw_unstand$neighbours[[Kerry]] <- red_lw_unstand$neighbours[[Kerry]][-Clare_in_Kerry]
red_lw_unstand$weights[[Clare]] <- red_lw_unstand$weights[[Clare]][-Kerry_in_Clare]
red_lw_unstand$weights[[Kerry]] <- red_lw_unstand$weights[[Kerry]][-Clare_in_Kerry]
-summary(red_lw_unstand)
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 25
@@ -565,7 +565,7 @@
Measures of spatial autocorrelation
Prop_stdR <- lapply(vars, function(x) moran.test(eire_ge1[[x]], listw=lw_std, randomisation=TRUE))
})
+## 0.118 0.000 0.118
## user system elapsed
-## 0.116 0.000 0.116
res <- sapply(c("MoranN", "MoranR", "GearyN", "GearyR", "Prop_unstdN", "Prop_unstdR", "Prop_stdN", "Prop_stdR"), function(x) sapply(get(x), "[[", "statistic"))
rownames(res) <- vars
@@ -608,7 +608,7 @@
Measures of spatial autocorrelation
I <- sapply(Prop_stdN, function(x) x$estimate[1])[raw_data]
EI <- sapply(Prop_stdN, function(x) x$estimate[2])[raw_data]
res <- (I - EI)/wrong_N_sqVI
-names(res) <- vars[raw_data]
+names(res) <- vars[raw_data]
print(formatC(res, format="f", digits=4), quote=FALSE)
## pagval2_10 pagval10_50 pagval50p cowspacre ocattlepacre pigspacre
## 3.8836 1.4957 4.8276 4.6744 1.9003 3.1550
@@ -621,7 +621,7 @@
Measures of spatial autocorrelation
calculated above:
res <- lapply(c("MoranR", "GearyR", "Prop_unstdR", "Prop_stdR"), function(x) sapply(get(x), function(y) c(y$estimate[1], sqrt(y$estimate[3]))))
-res <- t(do.call("rbind", res))
+res <- t(do.call("rbind", res))
colnames(res) <- c("I", "sigma_I", "C", "sigma_C", "unstd_r", "sigma_r", "std_r", "sigma_r")
rownames(res) <- vars
print(formatC(res, format="f", digits=4), quote=FALSE)
oMoranf <- function(x, nb) {
- z <- scale(x, scale=FALSE)
+ z <- scale(x, scale=FALSE)
n <- length(z)
glist <- lapply(1:n, function(i) {ii <- nb[[i]]; ifelse(ii > i, 1, 0)})
lw <- nb2listw(nb, glist=glist, style="B")
@@ -710,7 +710,7 @@ Simulating measures of s
}
f_bpara <- function(x, nsim, listw) {
boot(x, statistic=MoranI.boot, R=nsim, sim="parametric", ran.gen=Nsim,
- mle=list(mean=mean(x), sd=sd(x)), listw=listw, n=length(x),
+ mle=list(mean=mean(x), sd=sd(x)), listw=listw, n=length(x),
S0=Szero(listw))
}
nsim <- 4999
@@ -737,8 +737,8 @@ Simulating measures of s
Prop_stdRb <- lapply(vars, function(x) f_bperm(x=eire_ge1[[x]], nsim=nsim, listw=lw_std))
})
-res <- lapply(c("MoranNb", "MoranRb", "Prop_unstdNb", "Prop_unstdRb", "Prop_stdNb", "Prop_stdRb"), function(x) sapply(get(x), function(y) (y$t0 - mean(y$t))/sd(y$t)))
-res <- t(do.call("rbind", res))
+res <- lapply(c("MoranNb", "MoranRb", "Prop_unstdNb", "Prop_unstdRb", "Prop_stdNb", "Prop_stdRb"), function(x) sapply(get(x), function(y) (y$t0 - mean(y$t))/sd(y$t)))
+res <- t(do.call("rbind", res))
colnames(res) <- c("MoranNb", "MoranRb", "Prop_unstdNb", "Prop_unstdRb", "Prop_stdNb", "Prop_stdRb")
rownames(res) <- vars
We collate the results to compare with the analytical standard @@ -786,7 +786,7 @@
## user system elapsed
-## 0.064 0.000 0.065
+## 0.066 0.000 0.067
res <- sapply(c("MoranSad", "Prop_unstdSad", "Prop_stdSad"), function(x) sapply(get(x), "[[", "statistic"))
rownames(res) <- vars
## user system elapsed
-## 0.081 0.000 0.082
+## 0.084 0.000 0.084
res <- sapply(c("MoranEx", "Prop_unstdEx", "Prop_stdEx"), function(x) sapply(get(x), "[[", "statistic"))
rownames(res) <- vars
-vars_scaled <- lapply(vars, function(x) scale(eire_ge1[[x]], scale=FALSE))
+vars_scaled <- lapply(vars, function(x) scale(eire_ge1[[x]], scale=FALSE))
nb_W <- nb2listw(lw_unstand$neighbours, style="W")
pre <- spatialreg:::preAple(0, listw=nb_W)
MoranAPLE <- sapply(vars_scaled, function(x) spatialreg:::inAple(x, pre))
diff --git a/docs/articles/nb.html b/docs/articles/nb.html
index f4fc7ea1..791ddd08 100644
--- a/docs/articles/nb.html
+++ b/docs/articles/nb.html
@@ -133,7 +133,7 @@ Introductionlibrary(spdep)
## Loading required package: spData
## Loading required package: sf
-## Linking to GEOS 3.12.1, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUE
+## Linking to GEOS 3.12.1, GDAL 3.8.3, PROJ 9.3.1; sf_use_s2() is TRUE
NY8 <- as(sf::st_read(system.file("shapes/NY8_utm18.shp", package="spData")), "Spatial")
## Reading layer `NY8_utm18' from data source
@@ -152,8 +152,8 @@ Introduction
Syracuse <- NY8[NY8$AREANAME == "Syracuse city",]
-Sy0_nb <- subset(NY_nb, NY8$AREANAME == "Syracuse city")
-summary(Sy0_nb)
+Sy0_nb <- subset(NY_nb, NY8$AREANAME == "Syracuse city")
+summary(Sy0_nb)
## Neighbour list object:
## Number of regions: 63
## Number of nonzero links: 346
@@ -185,7 +185,7 @@ Creating Contiguity Neighbours## [1] "sp"
+isTRUE(all.equal(Sy0_nb, Sy1_nb, check.attributes=FALSE))
## [1] TRUE
As we can see, creating the contiguity neighbours from the
Syracuse
object reproduces the neighbours from (Waller and Gotway 2004). Careful examination of
@@ -204,21 +204,21 @@
Sy2_nb <- poly2nb(Syracuse, queen=FALSE)
-isTRUE(all.equal(Sy0_nb, Sy2_nb, check.attributes=FALSE))
+isTRUE(all.equal(Sy0_nb, Sy2_nb, check.attributes=FALSE))
## [1] FALSE
oopar <- par(mfrow=c(1,2), mar=c(3,3,1,1)+0.1)
-plot(Syracuse, border="grey60")
-plot(Sy0_nb, coordinates(Syracuse), add=TRUE, pch=19, cex=0.6)
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="a)", cex=0.8)
-plot(Syracuse, border="grey60")
-plot(Sy0_nb, coordinates(Syracuse), add=TRUE, pch=19, cex=0.6)
-plot(diffnb(Sy0_nb, Sy2_nb, verbose=FALSE), coordinates(Syracuse),
+plot(Syracuse, border="grey60")
+plot(Sy0_nb, coordinates(Syracuse), add=TRUE, pch=19, cex=0.6)
+text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="a)", cex=0.8)
+plot(Syracuse, border="grey60")
+plot(Sy0_nb, coordinates(Syracuse), add=TRUE, pch=19, cex=0.6)
+plot(diffnb(Sy0_nb, Sy2_nb, verbose=FALSE), coordinates(Syracuse),
add=TRUE, pch=".", cex=0.6, lwd=2)
## Warning in diffnb(Sy0_nb, Sy2_nb, verbose = FALSE): region.id differ; using ids
## of first list
+text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="b)", cex=0.8)
par(oopar)
Similar approaches may also be used to read ArcGIS coverage data by tallying the left neighbour and right neighbour arc indices with the polygons in the data set, using either RArcInfo or @@ -302,20 +302,20 @@
oopar <- par(mfrow=c(2,2), mar=c(1,1,1,1)+0.1)
-plot(Syracuse, border="grey60")
-plot(Sy4_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="a)", cex=0.8)
-plot(Syracuse, border="grey60")
+plot(Syracuse, border="grey60")
+plot(Sy4_nb, coords, add=TRUE, pch=".")
+text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="a)", cex=0.8)
+plot(Syracuse, border="grey60")
if (!is.null(Sy5_nb)) {
- plot(Sy5_nb, coords, add=TRUE, pch=".")
- text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="b)", cex=0.8)
+ plot(Sy5_nb, coords, add=TRUE, pch=".")
+ text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="b)", cex=0.8)
}
-plot(Syracuse, border="grey60")
-plot(Sy6_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="c)", cex=0.8)
-plot(Syracuse, border="grey60")
-plot(Sy7_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="d)", cex=0.8)
par(oopar)
oopar <- par(mfrow=c(1,3), mar=c(1,1,1,1)+0.1)
-plot(Syracuse, border="grey60")
-plot(Sy8_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="a)", cex=0.8)
-plot(Syracuse, border="grey60")
-plot(Sy9_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="b)", cex=0.8)
-plot(Syracuse, border="grey60")
-plot(Sy10_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="c)", cex=0.8)
par(oopar)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 395.7 587.3 700.1 760.4 906.1 1544.6
@@ -429,15 +429,15 @@Distance-Based Neighbours## 4 1 1
oopar <- par(mfrow=c(1,3), mar=c(1,1,1,1)+0.1)
-plot(Syracuse, border="grey60")
-plot(Sy11_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="a)", cex=0.8)
-plot(Syracuse, border="grey60")
-plot(Sy12_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="b)", cex=0.8)
-plot(Syracuse, border="grey60")
-plot(Sy13_nb, coords, add=TRUE, pch=".")
-text(bbox(Syracuse)[1,1], bbox(Syracuse)[2,2], labels="c)", cex=0.8)
par(oopar)
Distance-based neighbours: frequencies of numbers of neighbours by @@ -476,7 +476,7 @@
dsts0 <- unlist(nbdists(NY_nb, coordinates(NY8)))
-summary(dsts0)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 82.7 1505.0 3378.7 5865.8 8954.3 38438.1
If the areal entities are approximately regularly spaced, using @@ -506,14 +506,14 @@
Sy0_nb_lags <- nblag(Sy0_nb, maxlag=9)
-names(Sy0_nb_lags) <- c("first", "second", "third", "fourth", "fifth", "sixth", "seventh", "eighth", "ninth")
+names(Sy0_nb_lags) <- c("first", "second", "third", "fourth", "fifth", "sixth", "seventh", "eighth", "ninth")
res <- sapply(Sy0_nb_lags, function(x) table(card(x)))
mx <- max(unlist(sapply(Sy0_nb_lags, function(x) card(x))))
nn <- length(Sy0_nb_lags)
res1 <- matrix(0, ncol=(mx+1), nrow=nn)
-rownames(res1) <- names(res)
+rownames(res1) <- names(res)
colnames(res1) <- as.character(0:mx)
-for (i in 1:nn) res1[i, names(res[[i]])] <- res[[i]]
+for (i in 1:nn) res1[i, names(res[[i]])] <- res[[i]]
res1
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
## first 0 1 1 5 9 14 17 9 6 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
diff --git a/docs/articles/nb_sf.html b/docs/articles/nb_sf.html
index 865e1ecd..736105d0 100644
--- a/docs/articles/nb_sf.html
+++ b/docs/articles/nb_sf.html
@@ -170,8 +170,8 @@ Data setsf_bna$AREAKEY <- gsub("\\.", "", sf_bna$Primary.ID)
data(NY_data, package="spData")
key <- as.character(nydata$AREAKEY)
-sf_bna1 <- sf_bna[match(key, sf_bna$AREAKEY), c("AREAKEY")]
-sf_bna2 <- merge(sf_bna1, nydata, by="AREAKEY")
+sf_bna1 <- sf_bna[match(key, sf_bna$AREAKEY), c("AREAKEY")]
+sf_bna2 <- merge(sf_bna1, nydata, by="AREAKEY")
sf_bna2_utm18 <- st_transform(sf_bna2, "+proj=utm +zone=18 +datum=NAD27")
st_write(sf_bna2_utm18, "NY8_bna_utm18.gpkg")
@@ -229,7 +229,7 @@ ## GEOS GDAL proj.4 GDAL_with_GEOS USE_PROJ_H
-## "3.12.1" "3.8.2" "9.3.1" "true" "true"
+## "3.12.1" "3.8.3" "9.3.1" "true" "true"
## PROJ
## "9.3.1"
Let us read the GPKG file with valid geometries in to ‘sf’ and ‘sp’ @@ -258,7 +258,7 @@
## user system elapsed
-## 0.0934 0.0048 0.0995
+## 0.0989 0.0056 0.1050
Using spatial indices to check intersection of polygons is much faster than the legacy method in poly2nb. From spdep 1.1-7, use is made of GEOS through sf to find candidate @@ -306,7 +306,7 @@
try(NY8_sf_old_1_nb <- poly2nb(NY8_sf_old), silent = TRUE)
-all.equal(NY8_sf_old_1_nb, NY8_sf_1_nb, check.attributes=FALSE)
## [1] "Component 57: Numeric: lengths (4, 5) differ"
## [2] "Component 58: Numeric: lengths (5, 6) differ"
## [3] "Component 66: Numeric: lengths (7, 11) differ"
@@ -347,7 +347,7 @@ Contiguity neighbours from
geometries in the same ways as before imposing validity:
try(NY8_sf_old_1_nb_val <- poly2nb(NY8_sf_old_val), silent = TRUE)
-all.equal(NY8_sf_old_1_nb_val, NY8_sf_1_nb, check.attributes=FALSE)
+all.equal(NY8_sf_old_1_nb_val, NY8_sf_1_nb, check.attributes=FALSE)
## [1] "Component 57: Numeric: lengths (4, 5) differ"
## [2] "Component 58: Numeric: lengths (5, 6) differ"
## [3] "Component 66: Numeric: lengths (7, 11) differ"
@@ -356,7 +356,7 @@ Contiguity neighbours from
The neighbour sets are the same for the old boundaries with or
without imposing validity:
-all.equal(NY8_sf_old_1_nb_val, NY8_sf_old_1_nb, check.attributes=FALSE)
+all.equal(NY8_sf_old_1_nb_val, NY8_sf_old_1_nb, check.attributes=FALSE)
## [1] TRUE
@@ -383,7 +383,7 @@
if (unname(sf_extSoftVersion()["GEOS"] >= "3.9.0"))
- NY8_cic_sf <- st_cast(st_inscribed_circle(st_geometry(NY8_sf), nQuadSegs=0), "POINT")[(1:(2*nrow(NY8_sf)) %% 2) != 0]
We need to check whether coordinates are planar or not:
st_is_longlat(NY8_ct_sf)
system.time(for (i in 1:reps) NY88_nb_sf <- knn2nb(knearneigh(NY8_ct_sf, k=1)))/reps
## user system elapsed
-## 0.0182 0.0009 0.0193
+## 0.0181 0.0008 0.0191
Legacy code may be used omitting the kd-tree:
system.time(for (i in 1:reps) NY89_nb_sf <- knn2nb(knearneigh(NY8_ct_sf, k=1, use_kd_tree=FALSE)))/reps
## user system elapsed
-## 0.0174 0.0011 0.0186
+## 0.0180 0.0013 0.0195
dsts <- unlist(nbdists(NY88_nb_sf, NY8_ct_sf))
-summary(dsts)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 82.85 912.85 1801.11 3441.04 4461.26 17033.11
@@ -443,7 +443,7 @@Distance neighbours
system.time(for (i in 1:reps) NY810_nb <- dnearneigh(NY8_ct_sf, d1=0, d2=0.75*max_1nn))/reps
## user system elapsed
-## 0.0310 0.0007 0.0319
+## 0.0324 0.0007 0.0334
By default, the function uses dbscan::frNN()
to build a
kd-tree in 2D or 3D which is then used to find distance neighbours. For
small n, the argument use_kd_tree=FALSE
may speed up
@@ -453,7 +453,7 @@
system.time(for (i in 1:reps) NY811_nb <- dnearneigh(NY8_ct_sf, d1=0, d2=0.75*max_1nn, use_kd_tree=FALSE))/reps
## user system elapsed
-## 0.0166 0.0009 0.0175
+## 0.0176 0.0008 0.0185
## user system elapsed
-## 0.0245 0.0000 0.0246
+## 0.0251 0.0000 0.0252
For this smaller data set, the legacy approach without spatial indexing is adequate, but slows down as the number of observations increases:
@@ -493,12 +493,12 @@
system.time(for (i in 1:reps) pts_ll2_nb <- knn2nb(knearneigh(pts_ll, k=6)))/reps
## user system elapsed
-## 0.0186 0.0001 0.0188
+## 0.0202 0.0000 0.0203
The WGS84 ellipsoid Great Circle distances differ a very little from the s2 spherical distances, yielding output that here diverges for two tract centroids:
-all.equal(pts_ll1_nb, pts_ll2_nb, check.attributes=FALSE)
all.equal(pts_ll1_nb, pts_ll2_nb, check.attributes=FALSE)
## [1] "Component 52: Mean relative difference: 1.466667"
## [2] "Component 124: Mean relative difference: 0.0251046"
@@ -546,7 +546,7 @@Distance neighbours=0.75*max_1nn_ll))/(reps/5) }
## user system elapsed
-## 0.0430 0.0000 0.0435
+## 0.0445 0.0000 0.0450
Alternatively, spherical distances can be used with
dwithin=FALSE
and s2::s2_closest_edges()
;
although running in similar time, s2::s2_closest_edges()
@@ -555,9 +555,9 @@
system.time(for (i in 1:(reps/5)) pts_ll5_nb <- dnearneigh(pts_ll, d1=0, d2=0.75*max_1nn_ll, dwithin=FALSE))/(reps/5)
## user system elapsed
-## 0.0265 0.0005 0.0270
+## 0.029 0.000 0.029
-if (packageVersion("s2") > "1.0.7") all.equal(pts_ll3_nb, pts_ll5_nb, check.attributes=FALSE)
if (packageVersion("s2") > "1.0.7") all.equal(pts_ll3_nb, pts_ll5_nb, check.attributes=FALSE)
## [1] TRUE
Using s2::s2_closest_edges()
respects
d1 > 0
without requiring a second pass in R, so is
@@ -568,7 +568,7 @@
## user system elapsed
-## 0.027 0.000 0.027
+## 0.0275 0.0000 0.0280
Using s2::s2_dwithin_matrix()
requires a second pass,
one for the lower bound, another for the upper bound, and a set
difference operation to find neighbours in the distance band:
## user system elapsed
-## 0.0720 0.0005 0.0730
+## 0.0740 0.0000 0.0745
-if (packageVersion("s2") > "1.0.7") all.equal(pts_ll3a_nb, pts_ll5a_nb, check.attributes=FALSE)
if (packageVersion("s2") > "1.0.7") all.equal(pts_ll3a_nb, pts_ll5a_nb, check.attributes=FALSE)
## [1] TRUE
Setting use_s2=FALSE
falls back to the legacy version,
which uses symmetry to reduce time:
system.time(for (i in 1:reps) pts_ll6_nb <- dnearneigh(pts_ll, d1=0, d2=0.75*max_1nn_ll, use_s2=FALSE))/reps
## user system elapsed
-## 0.0094 0.0000 0.0095
+## 0.0101 0.0000 0.0102
Minor differences may occur between the legacy ellipsoid and s2 spherical approaches:
-all.equal(pts_ll5_nb, pts_ll6_nb, check.attributes=FALSE)
all.equal(pts_ll5_nb, pts_ll6_nb, check.attributes=FALSE)
## [1] "Component 20: Numeric: lengths (6, 5) differ"
## [2] "Component 28: Numeric: lengths (7, 6) differ"
## [3] "Component 112: Numeric: lengths (109, 108) differ"
@@ -617,9 +617,9 @@ Distance neighbours
system.time(for (i in 1:reps) pts_ll6a_nb <- dnearneigh(pts_ll, d1=5, d2=0.75*max_1nn_ll, use_s2=FALSE))/reps
## user system elapsed
-## 0.0096 0.0000 0.0096
+## 0.0099 0.0001 0.0100
-if (packageVersion("s2") > "1.0.7") all.equal(pts_ll5a_nb, pts_ll6a_nb, check.attributes=FALSE)
if (packageVersion("s2") > "1.0.7") all.equal(pts_ll5a_nb, pts_ll6a_nb, check.attributes=FALSE)
## [1] "Component 20: Numeric: lengths (6, 5) differ"
## [2] "Component 28: Numeric: lengths (7, 6) differ"
## [3] "Component 112: Numeric: lengths (62, 61) differ"
@@ -668,9 +668,9 @@ Contiguity neighbours for spherical polygon support
system.time(for (i in 1:reps) NY8_sf_1_nb_ll <- poly2nb(NY8_sf_ll, queen=TRUE, snap=eps))/reps
## user system elapsed
-## 0.1571 0.0022 0.1613
+## 0.1574 0.0024 0.1604
-all.equal(NY8_sf_1_nb, NY8_sf_1_nb_ll, check.attributes=FALSE)
all.equal(NY8_sf_1_nb, NY8_sf_1_nb_ll, check.attributes=FALSE)
## [1] TRUE
diff --git a/docs/articles/sids.html b/docs/articles/sids.html
index 837ad497..39dafff8 100644
--- a/docs/articles/sids.html
+++ b/docs/articles/sids.html
@@ -125,7 +125,7 @@ We will be using the spdep and
spreg packages, here version: spdep, version 1.3-2,
-2024-01-02, the sf package and the
+2024-02-06, the sf package and the
tmap package. The data from the sources referred to
above is documented in the help page for the nc.sids
data
set in spData. The actual data, included in a shapefile
@@ -151,21 +151,21 @@
sf_use_s2(FALSE)
-plot(st_geometry(nc), axes=TRUE)
-text(st_coordinates(st_centroid(st_geometry(nc), of_largest_polygon=TRUE)), label=nc$FIPSNO, cex=0.5)
We can examine the names of the columns of the data frame to see what it contains — in fact some of the same columns that we will be examining below, and some others which will be useful in cleaning the data set.
-names(nc)
names(nc)
## [1] "CNTY_ID" "AREA" "PERIMETER" "CNTY_" "NAME" "FIPS"
## [7] "FIPSNO" "CRESS_ID" "BIR74" "SID74" "NWBIR74" "BIR79"
## [13] "SID79" "NWBIR79" "east" "north" "x" "y"
## [19] "lon" "lat" "L_id" "M_id" "geometry"
-summary(nc)
summary(nc)
## CNTY_ID AREA PERIMETER CNTY_
## Min. :1825 Min. :0.0420 Min. :0.999 Min. :1825
## 1st Qu.:1902 1st Qu.:0.0910 1st Qu.:1.324 1st Qu.:1902
@@ -324,7 +324,7 @@ Getting the data into R
-all.equal(as.data.frame(nc_sf)[,1:14], as.data.frame(sids_sf)[,1:14])
+all.equal(as.data.frame(nc_sf)[,1:14], as.data.frame(sids_sf)[,1:14])
## [1] "Names: 12 string mismatches"
## [2] "Component 4: Modes: numeric, character"
## [3] "Component 4: target is numeric, current is character"
@@ -342,7 +342,7 @@ Getting the data into R## [15] "Component 14: Attributes: < target is NULL, current is list >"
## [16] "Component 14: target is numeric, current is sfc_MULTIPOLYGON"
-all.equal(as.data.frame(nc_sf)[,1:14], as.data.frame(sids2_sf)[,1:14])
all.equal(as.data.frame(nc_sf)[,1:14], as.data.frame(sids2_sf)[,1:14])
## [1] "Names: 12 string mismatches"
## [2] "Component 4: Modes: numeric, character"
## [3] "Component 4: target is numeric, current is character"
@@ -360,7 +360,7 @@ Getting the data into RThe spData data set has some columns reordered and a
surprise:
-all.equal(as.data.frame(nc_sf)[,1:14], as.data.frame(nc)[,c(2,3,4,1,5:14)])
+all.equal(as.data.frame(nc_sf)[,1:14], as.data.frame(nc)[,c(2,3,4,1,5:14)])
## [1] "Component \"NWBIR74\": Mean relative difference: 0.04891304"
so a difference in NWBIR74
:
@@ -399,8 +399,8 @@Getting the data into R## 37055 37095 ## 3 disjoint connected subgraphs
-plot(st_geometry(nc), border="grey")
-plot(ncCC89, st_centroid(st_geometry(nc), of_largest_polygon), add=TRUE, col="blue")
plot(st_geometry(nc), border="grey")
+plot(ncCC89, st_centroid(st_geometry(nc), of_largest_polygon), add=TRUE, col="blue")
Printing the neighbour object shows that it is a neighbour list object, with a very sparse structure — if displayed as a matrix, only @@ -411,16 +411,16 @@
region.id
of
37001
can be retreived by matching against the
indices. More information can be obtained by using
-summary()
on an nb
+summary()
on an nb
object. Finally, we associate a vector of names with the neighbour list,
through the row.names
argument. The names
should be unique, as with data frame row names.
+ncCC89[[match("37001", r.id)]]
## [1] 11 26 29 30 48
-r.id[ncCC89[[match("37001", r.id)]]]
r.id[ncCC89[[match("37001", r.id)]]]
## [1] 37033 37081 37135 37063 37037
The neighbour list object records neighbours by their order in relation to the list itself, so the neighbours list for the county with @@ -510,7 +510,7 @@
-hist(nc$pmap, main="")
+hist(nc$pmap, main="")
One ad-hoc way to assess the impact of the possible failure of our assumption that the counts follow the Poisson distribution is to @@ -609,7 +609,7 @@
nc$ft.SID74 <- sqrt(1000)*(sqrt(nc$SID74/nc$BIR74) + sqrt((nc$SID74+1)/nc$BIR74))
-stem(round(nc$ft.SID74, 1), scale=2)
##
## The decimal point is at the |
##
@@ -652,9 +652,9 @@ Median polish smoothingmFT <- sqrt(1000)*(sqrt(mSID74/mBIR74) + sqrt((mSID74+1)/mBIR74))
# mFT1 <- t(matrix(mFT, 4, 4, byrow=TRUE))
# wrong assignment of 12 elements to a 4x4 matrix detected by CRAN test 2021-05-22
-rc <- do.call("rbind", lapply(strsplit(names(mFT), ":"), as.integer))
+rc <- do.call("rbind", lapply(strsplit(names(mFT), ":"), as.integer))
mFT1 <- matrix(as.numeric(NA), 4, 4)
-for (i in 1:nrow(rc)) mFT1[rc[i,1], rc[i,2]] <- mFT[i]
+for (i in 1:nrow(rc)) mFT1[rc[i,1], rc[i,2]] <- mFT[i]
med <- medpolish(mFT1, na.rm=TRUE, trace.iter=FALSE)
med
##
@@ -681,8 +681,8 @@ Median polish smoothing
-mL_id <- model.matrix(~ as.factor(nc$L_id) -1)
-mM_id <- model.matrix(~ as.factor(nc$M_id) -1)
+mL_id <- model.matrix(~ as.factor(nc$L_id) -1)
+mM_id <- model.matrix(~ as.factor(nc$M_id) -1)
nc$pred <- c(med$overall + mL_id %*% med$row + mM_id %*% med$col)
nc$mp_resid <- nc$ft.SID74 - nc$pred
diff --git a/docs/articles/sids_files/figure-html/unnamed-chunk-3-1.png b/docs/articles/sids_files/figure-html/unnamed-chunk-3-1.png
index a1e842960b7a86396deb6b4ffd4fe05df2c0ada2..c0569c25afb63e0a30f2b44d492743f757556a59 100644
GIT binary patch
literal 116768
zcmeFYXH*kw+cr!G6%Y)HQlt~%M_Q<%Y)V8yiXfpW(tGbE2uKON
zw-BUC?*s@TZ`}93pYPxI`(5jOSg>Z8nQLaw<2tY7tONTCtU_~>=_U~o5sm6AkPZ_n;{
zMLkIR4$Ak7pnmH4;r_s2$}I~RCnt%Ti89$;Mp5?VI|ClX^mmmRjqi5yH@$!MH7r|z
zzaWHoRG*>aZ({1#_ZUX>8I(V3jTii}ttc#v%e4lKXux`ahz-+{;2#A~L@3_2$48
z{yA@hV19|Z@gIu(_)&q#@_{#y_cO5tjcNaa9I2Gu-?&fz3>~wVW1M-Q;y?VJ*y3aV
zGtg!0f9|-K>BB+N0*~0#_I>^@&rUFtw*0#*{*C%T?LYfte~Tdm{-t6-KT-L=bD{a@
z|Bj0bF23&}SbT54?g*J{bA7)1?_0y^C?J&EUa#YbYvrDP_2K!C0Gd?l;;;TktZ?6o
zMQZ1nHBTl4|B~<$c!~D4-4ckpSN%U%@IOY#7%60~9d8heP{#5t_)@wf1m{2hM}Be2
z`%A2r^CXt!bpi{%{0DcIrs7Zjdut<)Ui`T+Q48Rl)XyE`zf=iYFY7S>7sHSTEi3D9
z5Ss|w2&_CwOL9Tqf6>Kf5icHAFr2d=EdL0Cl4jmCdz|p_-%(5f$@ceKTSEjKX!jnH
z%Myh*|2tE2kHdYv_J{*-Li~TVp^3FfuVw#d5d4wt0WF&0z9D)YwO@IS+J
zGzn*>TSLY}R2cWL^am_UMgQOHcZr4jO67%iMw>MAH8kImB7#s;`2LIFDS~>jWB#TT
z>rACP`T^oNlv!nzX|LAgpxFq4ClJJ)(
z;{QG}C9&k+-Txmw?X0IZlHLOUJ3a6~3P{vjCJ6n3){?=J&yp4*`1?-Mf2K93&o;#O
zATMM*6cD11pbWW9%1mNM;>K(~{jb`q?#kY>C(0sQzq@Bi7AP686Uy~_x`0K!Wu5CK
zivB;Z@^QirCN=7*f{ItQrPyC$B`zVmCN*gtw62qpwuECSU&|Dsx4o;dB5aOv|GvTC4T5IxLte`jAs@4R#$7EOG;#e4
zVc8b~H!RQuzAxy+-gqs!WW5{&eRjaTB!igc+Rpe_LYNKVp$4S4$-%^#S*V1SkZHPV
z@_$#2jv@5{ZO<;^NeE4-ONc(mWzke+mhWGH69
z2)WR;T+_)QWXNf?sfdBLh*S%8L*#wdYW{7f$zg-I(#=ZfpMl~M)CaLyHi$FCON2RM
z4PgVh!>`Izwn03Om=1CCxmX}a$o*!aDr<6kR`K`DZCCMm#!36ZVqxLKv;vp-oX6%45dsc!p6v)&?!_Fo4v
z%#cfwLhRDk_K3|u)KR5XveT6p-XC<5NT&V&r`SrvgL28gGJ(<
z@fH8|QPbBA^F2+**Yw8dGI9gD*@2kqtmEjnVc!0Mt*s_XX3nnO&83|y`p6P)3pD(P
zi}%_1>#CgDFO-6ubs1}>b9us?FAFsWmnqyutgu1)F*YWp-vgmyDZyE`>(^y1XeVE^
zG@><>Il9kunvT(Xdg3Zb0x3*dM<&dDrs}P5yZS~GeGaGbVg@Ic9!fcM=
z3$ChBZ`W*+iJ3|8+-GZrxWXU%Ywf1E#l6ZL#bD>Zu7iZ!IkkUo+sqy_cZqhhg<9Tq
zg$^W(%rir7)3cc>x5N$e#ERy^-h}v8)N$Few-pDg0s}gAk%k%$!mZZa&gK>Cb*Jn
z4U|m^T~7>wP2Sl5#d)(TH|8(%{a4ulw?d@LR}TvK`7mhv&G@6hDae
zf6GzLaA9$#i1Fxfo12(cfYd{XsOYc3w;`U9)5@x2sNN<-brm{*Z}_p`CVQ|1@uJ5u
zd}qV0aEM8j7;YLh_5G7rVUa7+zX4ASI@l2AG}5Ql1Jm~uJj{=+hq8K8u76@U5VEbH
z(P<8G&u2g61tHFE7KPVuwsn&R(`4k_c~a^s&`Y}0x|kQbVZDuYRt_kgbMuCZ83v0?
zH}0Q26F}E9ozrV4Z{B=uX=|HnU{h3wZQWZDVILYNwou<5CUSYWeLAI{cAGes6J7uC
zm9;$MlSJfhJ$BY1Rk`cYT7V6vHSyKWS6REYciZa&O13Z)XdmM-q770Rk}U368a7Q9
zCAqj6+8vz+RMdTXT)~g+mJF&(5;KxOW+{>x0BH;S@k6Pqp}F;%LBDGk>x%=@n
zK;@zTT?z<~YWSfbCmn`sNsDVwT?W$|JoNc(YVIm9ai7^ng};Pz`0t&ldDwKo{U;)v
zIp0#dL-qr0DMGQ{P!X4$H=S6!Alxl?xyc9ZPNoa!O^LLPyVWn|X*zaK5Ah@tz?TKM
z+bEl2TK)hZj!A%oFl8N2ofD1KgqdXOk5*d--{bqaXbuu{NxqG@fq0
z_%LVOCMNe}Pceh>jq=onOs|T%1~EQC?)s@Mg@p{c?K~>dTrv=FM6|D{WQ0bIhJtL@
zgQ5{AJhekBo=jdTFnLiS$9w4K0|Z;hFo4kz%>5EzDTDBKrMY4|Z1*G8&a*bcXJvcB
z?)M3Iz1k2V(INwrWKy@ZR_y!~rwsJvIJ0X8;J0_9V|Muno*W*HfrZCt)fulZIH#^S
z)zpg&8h!wdR)X+x_q|oPRX+!VobxRoZ^e53wE%o*czAxVOL;;mI{btuq9^&2ox2|D
zy%Oz)4Te?iLqp-sN4_R`lb6Yi$$kN6CS#v4jn3qUCKzUK728eWWpCM_3DyTy>e17<
z3vJ(=JRlv;{MlA)&0|gHC-)L4H})l$PK7wDHuqJ@WF@%zd7YKWEPY8HQtVL3CKvd&
zi5zeol?>pbe{JL8LhUs4k7hqLP>ngt)O`t@mlN)64mf?j?ES%s?|b7i%MoscG_q!B
z8RLQBD6MA81Z9>6JSp8c-%RNSwa%SnJnVjU{iMnpDv1psK6rLCJdijq>znA@$XvZ4
z86ERf(yrZ2^owxu4?7s-!
zaip-<_PotfqvMCelCzclvIoeZZ#cSS)GOQ5(n`wW;w-A9G^1+m9^Kjj9rzQfd|T*ZWz~n7IOxcwmHnc?Im4gN?+?Dho9gY?m!sB6~Gsi8N+Km
zqZX^!+4<=GJ$4~7VIDl5<=N8dXP5L+Ol*k=T~8Q35!DW#QtUioz2>9OvZw3#Rc_nf
zUeSF#0i}REJ=-FKjNh%3Hl$!oaFe2J{q~T7;FF2BStvPw*Lvu!;t1y*W@(I(>W#Wz
zBw1dwlj>E**8&t3d?JUT;)Ix{`(08m{r=&(C3p;3Uj6o2s1LK>?Dc!a>{C@XuiO(}
zf{E@gf)~;fSjvEP9(*&!zaArsD-D4N(^^qS8{!KJ$24P%
zt2iO{89P^-y~W1lteFFosHKUeMe@Lk;cb6@KZ@Vn$>6z&Km~8c
z!$q30z$M>Rd?F#}VP&Rg)JNaC`f_zHf+~7B5hxN+Ql!^OX+I==t$T)vsmd}GJy#XE
zOMP$_D)V&|OM37MzI{+z1^$JZrx?#ro4>=ueOdWdJ+(j4!c#qAJ)7Y*OV6*4QZuC_
zN2xCbP*WRu(Sf}QiQOFy>xLZXj}ajAeswh@E|qG*6-4$yTFCRMDDvd9%TxwfdX7AF
zwxCPx%!{A7hD)>)c|3jbD{lz%~t1dSt5W+Crj>3)YONR^#
z^T{QR9IM7G;*2r0oZ1%YShuIxIx>5j>O8)7fEFnfMFxc$%g!QcHLu
zcQ?+1x*HxO6||!MtLyJK+>|H!Rv?@9-K
z_~b9?v2apLJFMPgwbn<+NZ$i>f>ntUDR3_f%G`PSg@dFsQLjm3EQIB`Gj$t?dZqFG
z)PqtHg8gu-vWa5`*+s2RcR8JJ707U^!rurGv}T&X@QeyvmlpbVz7IXbaM*jpa;d2R
z2SMcAe2Eu6Vdr|+NqAhXXn%^`Hm|+e%nR>(q-x1)-oJ!$BG}bzE1ESA0lpEneh>~-
zAT!qCx_T!BMp8OuZwvoR8g+NOofw!#w9^n8TxT*fCg8{^7!I}Ws_-<6cR^RrTh~PrIk`F#N
zcTH<}4+^|^6xoavoqLqpQ4AATb~Z`|bgHm{sE;?rH8=;iW{4h2_0rOQAw-ExmfQ66
zzu!#x?Wd%~?*lDO4$qTk->v~xKUQJug-{C0)~U%ro>_6zU>qE&XH2^Q7xD~9TU*aQ
zffn`WWQy}4OrP#}mYHN-1xWA>ab9AN2@5nF_=(c|rI2RBvw~Zq-9$$Tw8+@+S~{%v
z9ql*{?qnPMW*s0fuA;RX`jvNK)bhxQQ*;a*FNaMLt9T}gpt=Eh8{>|@FKw_nOX
zrPxBA|GIYIeC3
zS928}j^+m!KP7D$Bm0wxTXEvMAa-XOn-UfvyG$$9@l91Wwxs2VBnI2n$MiZG?V0)y
zaaH&9KixqqZ8sfQa)6WHXKEgo^m`}ssXytJynj_PwX!~V(hgoO$T6uY<;A4MJbd-=
z{0qYdU68z84eIyrO^`y*wIo#v;lqH8r7T6lzKDsEI^^Ie@-C>X_J|98t`yX
z@kFH7Xy5Ei@hGuRRnj5xV_R*_v`v)Yp}5oe8Q=?YsmG^j+*u{I-`<}r?yJF-^7A*`_1*W(45E!T#|ra)T-=t$3X5uDhm9TzwSnwa$?UP%v
z1nZKtFmlS{#~(x-EREajr8FGwUoe?Wv+ZtBrCT#EzpxeAE}jR=IFJpTp0o>ZX|G}r
zg`DK=3^+v{tQOu^db9!5vJp21_YYFESMI>h9B{9c78{}Ugg`2t*G~sf>E~6S{c8+
zX_&vrZjEEH8dIL1;zyev^P6+~5(M8>yJ*mVg+RFxdP$UgNJ|a*#>wXB#MOFI5hCGy{~vX7
zp(1oZDZ2ep!)OdP>HXSvoo;?9Ke8%o&3LGuU6{H=D?L3^<@&<-#nrHBBT+-G7g_?u
zz9b@)`6Q0Yz)j0c9Di>RFe1pw3vFVfUy#o_SI~2)*)$a!n4T_bMbo3pEv8aRvT!Zc
z-771Pc$91WHCH%4iW|v|S^Q1mahRtoG&1{Yeu%-{5yKjK4b{3wpXlbQtYLNbHe1*|
z{}e}Sl#UCw>0asd;ywT^WfB|}YHf8WkEl-0jwzz!>6{#|lgzGLGRX$XKE-DzuGIhb
z0?1D+Kl7e4%+CPpupbmX%t8`TDZ$*GiinS
z?IiiZo1xroqr`3P%xCvakWIb~(hP)GiyT1i$t3l*l&d{=s4oK}bKxJt(;t_x#qwPG
zbr;CI$n2=UAexTvCn>dRSR-S%n5OAcuIhv#eyhRTk-KdwD@m<*&gC%UIUq{4CCa8X
z7Qt8CT)JP==6a^kEit}rMkshSJ2boEP}b^2>Tu^IRJ+UOH6JI$qU}qQL*3Qoz(i4L
zQ;^z&Y61JKRl^0lpp+{ko*+X7HI@x;nZ`XRKymV^G4p0XTh++5W|!z>#wsJ<``txQ
zo8|RgdaCxnrxK6y_TAwlhx~9%#4#a=^ONw4JiKYuY;@YL4gJzl-29mbArKwj-V<9V
zP;K>xqmgp3@L7rqr6`*Xddl^7x9WjD#D*G1=&kJy@?NJcS4+>-ReDfr_(~`cK)FnF
zv@??JQ>|h+8zZIq3dkaPFIZ%8?-4FY9;Nxp>&;oqzV~5zwgs;bUvWT~xt&f*lM5~t
zg6f8S_GeXFKt5hpMmz;QuKWxJT<{cRr9<_SOB)H0L6FE>2{X(R_~P+_>z{UIvtTO
zFc&OWX%3zXii_0#k9Cl`9OD9tf48uHI!JuHCyzk
zDUFp*efIF3wa9+eTMqNp>5$>Eh)jxullRl5*M#hrCi|a6kA;Q0b4t6p=Y>Hk|Kn!o
z&lX}P2|c6I+oi69BMw%-?e!l^W(byZiz#L%%Gm+o;T=xVN4LDyUbl(_4ADGI*0hzW
zR5wNecEGz{<5LQRDN5?{=?JidhFRb}@-l(ic@J>E<0uP2_0q&wf>*&f6a;
zyBl{Z{mm|8zfZFIEar?u_?o-4O@%qZ?YO6~P6qURYg>H!uhttPTVa<;CxNE-qTu5qs&iwRaCMWOA%#
zeO$c9^WEsnD~r06faSfJDPdik4fgmSVyzlwX>A*YZoc$9^cKn%#R?y<`1VrD*@a&-
zV36H)S~lw*3mL1a=(Gf=)|~D=bEfON!^r=coIEq!GEJBn(B(pz$kAL^ZO08aPN%BZT4X)jr!_M7c#Ih
z$|WeDUZIr?q6_=tBrP-(ZnN5f^sO_hbrAYr-2Q21^)q|#j)Jxn5&6->uN>_3TWBfRR!8q&N#Z)l98)mVED8f$4+LAD3+t
zE0@w4vY2g9D_F*RK1}xOal#m`vlpD%KZUdAe4YwkhrX0=91=f#AZU3o6O#%o-n6b6
z+F9rvXbV=L7jNL-UMYvg(Cw1YcvsbMWZATtl{cq+_0AxA{<*26y~gZPL9XRshP_%T
zF#|@8A549)RN7Cg$-@+^``r>(L0-De#7HWHdt_W0ec$1b^R>beX4!4lKQ{kD=AEeM
zg9&!zjh6T4vJWY*GEx)cYA+w=Hl-&%9{pp`)v;;N?^t>%8x)!C{KxZ69z=lp)S?Pu
zPQH)#RlmxbvTgM1>dpi+c5Qv>pM$$dKv44L8ZAG;fG(doJGqbsB8!&C7FWWq&riF@
zPEsMh%FTscuP$RzS?96s;RN^{Wmr4&f_h2dB2CvxS5^Z|R`c+)2(-va5MH9(JIXng
zx!1{@Q3YB6{h!b+j$uZ&WV!OsOh9M(Gm`!8dwgw!OD!qDQ=e!vYBJ6>mrY%TwwXQ;T_310tv&9=2u;}M6BbKiubkQ4aYwvwpk_l+|R0Z?CQDP-#6VDsAE1H
zcL|XsKW0Kac`%R@%RFw@bAHoGW^{{0Y@)ndPWX%q8g$#f+6)zB;Bhj3d=pa-m-4OO
zx%O1sy$5V;{<{EXeih4zK?c}i^4VSHLLOpI)gamR*c?jz%|)n@~7&kaboGcE{SmwPs>2QwvE>
ztMsLXs<{JN%pW(}1CPQl#;dXS61C1Jzxx09WYH2XE7e+aVsw*QJVYeXE_{KcWV*A-
zKXoRzVEl%UWMb2VoUs;gIL@(p;A>SPquusmwHIE;mZEA~V>N{&YTo;|`fNp2c}(W>
zliDq-S13W3ULb4W)9TwM^*4!fuB?{M~u!h+bM2y=S2-(Rerh(JI$BsNFMx&43A_vYm=OZT>if
zs?tRH*Lr#>W&M4JdahuTT3sIVPp`^-+9axPYkX%XsHTFJgSsdyL|rp!BWVdw*coamZ~^&$qfQ<~;*|1l51os1xU~I9$Rb
z*6V*m92_P-Ya;!nz2yLd#5QDJfNVX&rBoj^@ii-GveyDGV5_+`v{I*@Rw>y&z!_A{
zcwuRFQ>`3je1@!Y3tIoFCo{-aTREAyap0zp{(y}<
zA4wiD-1xDsLwwOb1fGQX3ywjS@~)vW0@KcZ{E&7PS(x6d-Alw}2#Z~XU*4c78YU|L
zX;B0?Aw;(tiNj3B8$Q2@urk#DsP#gMEEfk&Gnles|MitE<$jcc=#v8TGV=M`h4l$BCTW1F0Pm;Y|K^4sK9XFdWf0QF?2Tt=iba(r(R={XqAS3lB>o_w@m
z2zT(8-ru+hoA4dWT>6`lG$S9
zRG(%dy~lJTq0L_p_g2x^W2X}6#28gN1541*r>W_3dWF96KFlnisV7~cbY6^lcfzga
z_^kFkc?L|kKy*y0h}ApC>ctdUbv)EbK7UBwNrrP?x;6-q0UXlYy$*Q=k%3>yx%$Rm
z2y}_6d9*LM>MRKN&CS->{cg6|mzt
zL+ZevK;PFc8t4n+Gk>hqi&e|Jh!rGmyfIF>^3JaKGp*n3c~n)NrDDxAgq*^#b5bXY>S?AyD)QBDF|@`-&_ZLoN50rd#P8{WkyZ0$V76>m
zjM}b+=9R*}mh}{y0p!{wdsWG=J<#B6rqe-=&+SmAC1fTt*JqeQJ{^}?7qn}!J*~cP
z`s`83>=Qzv(m)pG;M*;40|<&ao7F+W5_c&?ldN_A$@NumT&SuQ-EB@PpfMfzt1;wE
z$w0;CtreA95iA+ZeN7w~CB-q5U2&CWPg_GPoyb>OE%gTZ6d_4}A}!SmE|3
zy0~E~_I_v89QIt(1a0U*%P^JED-xttR{#-w(7i2_y3Uj^O`p28=F>blmN&9qy{R?g
zAFwR9n7WmC($o=YfXscd-;NT{0eLOBc9D>9e!wjc(^)LJj*X^MZ!IU-fWflC8qg}$TM
z`nbBkjk$@$IAwWrJZYa^n|DEoK23HaQm|SNzA&!IN~*tNZOaVIu3;Px^k^M-3UY90
z9M^7turz6;4)Dnl7}HO_o9y8G?vrd`MPRZ2%ax&j6y0wXSY-L#uf2*!^oFWv__AJ6
zxLK!!X8Ji}_=r>LGiNcB!HC(k3wdV-v7ZvBX)ek>)7q00n?Ez{3Or!4<146ZTCUW#
z28!gK9_nSluEy+KkD0ca&bHM7ch5Y7x>u%xDN3_mZ2t>tP{dWUp50M?u#GmO{%89F
zf=AeqtTkm{jqbM*(6M}&_z#w2)`=n<_yHf&upGnjE_}1h7tc98G4`D5Ee&$g
z5x;CszyNO>Oifw};3J`evASgM;Pk$=(J0!szxJNi%`IF@LR2j6Gz
z0&|+@mU-%z6B--&o(A{-Q-Yt#ls&(=VY?TTsp<*)Dgz4{e6*{+m}MIMp+Yk{{1iyk
z1E9QI(|}X>$t`$btd3xF(v$boxwUu{`qFOcgw)Hb&Y0;BA@=MDjQn?fe{XSf+0i>b
zDLvqQ8XPiPxaj7AflZV}7VbKp$%Xa>!P?#I-|1T-xOB!-YBfx>e&zauq}tB!8>2vMG|`Vkt8<(tRa+xAFI|hpb08^q4`2!@e$wRe4Zjs2BZ;
zGlfjNOS3cTw$1ip*rUOBREIb)U?{JzKba`b(bKfvs(EtiHki7U%1fWzMWpBU|mM
z6s71}9cws5fPD~sJ7oKe^CdsIws&Q@+B1QRU*=_6y9Cvo!9FFJ55e2$DI
zjpZ4_jCa4M5Un3tsXpC%VZpQ{c>>LS^x*k}iI*rPK0S0$#ALs?v#KN-ZyP4}-z%t47ZP+Mt`+)13%V}ct+o;K|Zd|(^
zG6UR3W2~C=(_iM4_a>N8Xg~Y-?+m)!n)myM@lP|qUlNpUrxkYo={jQTqXsO1wAylIA(+|7)hALOFDWoPft@48a
zJ+sFJgbvfWxIidqn=MbA;XMEQ-m9#ZScu$oToYH2*k-B+N)52FKuLc_CllSOdUgzBrQ0|0Jc)WjvA9K_xUAFTU
z)jne<_hd#qGq$YfRO=Lvz5~GVx7RQaawWU+4o?Zq$ES)$J?Rq>Mk&h!m1Z+Rd>JQK
z-o5Ij=^Ce)+J?4^)zpFR)w<*8eSI-yXTbXhxD3I(Mh1`{=(!BKHyIlA(PrFE7}Ykf
z%r==y_@bayJpDVNR976BmD>5j#E;kSdBc!?p~oNB_R>w>V=N)^ba>rIBg(tArDDhM
z@=AJ3io+D6&t-(EU>#VETi{dch;QOb`DuQMP6FT|@5c1+^m`
zXO#805PMi5BYo)XakTWdx|$0vS$0IHB`A#Q1uimB>cqP%IM*@ArWI^|MsPKuQ}#|1
z!t9rqt%Ry&s2r@V%toy0V{uNI5^7IkDE5{MY+GXQhRC+bs2q8avBBuUUqqGFr2oRt1w*p$@yJ;eCHQc
zT8%Yr-O#sCF$>gq(C01kV<%Av`YdBZ3a(2$k}KaK>kxW+BpG&L}Ow{4Y;f
z_*oWx=%=6pN!W(Ut!-u>+mW}uS|&j6AET~jDlLn{9|(GPaN{Um(^K)1cH+U@Gscz(
zNb7WL0~<;z^Y^h0K|ybIQ`hVD9CKPRqrlJ7XqRPI>Q^@ax`+qjc0G4<{kyRSHQt6(
zDfu-7)iKeuPj>Vy|Drrk=&eu@AvC8s&V1bBIP}i+C82+gBmS<=KNlhkF-||J?O5XV
zCWG9k=eE!c_<5*+nrMRE%
z=@Z8q{4+SV&F*%P9-$)66o4VxQ$GGM1LYQbCjkUgNpP2|VxDE6{D_R0|JeoeyScZh
zI1!VoM+gN975a^vpTjkBqt)!Dll_(4iVXlhQYUWBGC?Q<;0?e!n-KyR}IISae0
zZfOZ_N2mN-WKOWY*|B|nM(TTazUN453ZHo|JwFQZ0Hy&H*OD3o7%s4%)YQ#P{D)>i
zb?yZBs%7o8HA2cDNTyOmwS1*7Vq>R!=cVx7ZF)S#OxD>AVh8KFXHVzkL)I}+9c0!A
zu0L39ZY{98CjR+?XK=BdjyS|>YOg{|Dj)W%orZJXrk
zK`~rXn9+*tx79x=w6s`;Q!TcJzY@QEbO0SFF@=-F8q%Xvz0kKJQ6AYPYT=_vVN6My
zupqy1NKm1yDk((L@3|lJXtnG_fe)W$Ktu@op#eUY8KZZb@AQ1I<5H;!`X|WN^YFMH4pFesoH%o#n{8-)b-EnGC@9PoHn0XCq{AL0
z*_aPm2%BaFF1J&Nl7Mk<2$A_P(Z0qX1@S$GHM#ORVcv0$ztKyMpXO?gS-kzT+1s&h
zsGZg!=)2dmLz6hdO1NC=>2gZ_`|`7Q&f6-&!s`8$uG+E)p3Pn1w&|;ddpWS=P}ehj
zkotU(dLyq7DbPYr$V?RE=sE~D*24YKts*1y-SvK;iuUuC
zENd4OIPe~r_vM3n-Fi@`Eaz*>bdlw^t2j5h@4^yxKG>;lQ0X~^b{N;^*YTC
z#?&1C?!WPLPOC?ZujALA*n7K?KgMwvOzpzII7Sx%h#sSH2)FQVzHx53d~+H;_VypW
zIbE*={TzH$?NjWmD4&j@hy?6N$H;Oylr_JXZMGZBeciG7#}MZ1TNNAkY{1z>{fNuF
zIOT_`GxDkHy-J`=elN|9i=CT3R%0a}p%fIFp((<1VX{k?{SwF%U4EVshlF8mc+O?B
zO~9Sg)N-Ed`sq82&a~MRL&6p)lpcD(o-5WQcOoqZypz>G)Fo-pTIyR?Q<
z2`YRH3~OdUTu=Wvnoq!MP3@1I#f#(1`*S`1d<+MFat~k_yr;2
zH@3KeOjAMg1I43ZY`Ig7wj@GESzdB9CDnN4++mebM>xAF{e`*Kz+7WR&9zE=oU8!l
z9bJcLAD(yYv6?Z5skeu-v-%0nX?mTIhr8hdBrtQG_3EG20~D6(qV4tFCsD~C!^3;A
zGv&^W;vGhZ$&sSy9}unHaS_RQyRnZk4ojO|_zS@;OZRh`Llb$aufa5?mm0!5HA^h<
zj;>ldBQ^=Vd>V{>na@36npe@6R&VB2Z72}3TW~~4di+$Gko^$vPExC+DPKX=```Og
z(#e5IYV~AZj^@2{w}w@==fpX8Sq3HYF0(ILWeITdACT8R1gn!A35VNrJ&OPCQ+{`T8HcN+1N
zwy8AH8ovVQvO+IFgWC#sayj1T^~K3uo4v*96E7NUH-EU2+Vh)z%Qw<2dN`}W*c?pk
ztV1|=@2*9sPXu;94dn}bA`-ZhB6nsur={^sIaqd_NO-d{&e1l+&~uDrDMLJKz=BfN
zHYsC^z~c!M87lEvS9v;A8?iMirLEzCXm3i7u|pa7?Xqn+ufmOU;1gi@gZz*U_T1;l
zueKHCVw%<@<|<_O+cY!+j|8oPLr
zaPtE5+ndzFuspBt2|PeZq?+pq=}SA3ti6L4IT+5;e~%s%B5A&m&RI1cCoH=$es{Mm
zb6w-7#bqLiNxs+bW9Oi!ro*2}Erz9D0Wdha^OtAXcHh~bysPuEK8&U*iPO`X_4Cm>
zpY*X6T6dEJ45H)-iQ)Ou{fJp~{hc+MDLO^e-L8Nm_M?(jg(sdzy-i@Nb^CHXi(j|y
zBDE5=oaYO7e)QM}yI}Sz|oahKLI
zf^ZBnV#;3IE;%rAN+pxjp$F!K=c0z5Zp)kk$=OAIj|@5I4tH^zzuKC7#~G3pmO;~$
zRdQxZD7ZR2^f)ZgSB=fqOWw$XyG)y-8T?@HkE7}Pw?_gw90AM6UlLI7Pv33NNd{c+
z-`v*Z%E34trt3Yd)|FcUjs|K>xERw65}bU%+xG=p=d{_elL0UR=S|019nR{nDA+{(
z7iPoAfACbHpC5>U{j&XJSJhTgv2gDszs|CXhMuKBC%+CIi|9bWR#UZTIE~E9=#713
z2eIvh<`ISIT#odBX|rI;BmZQ{aoK>&j|g-K)VhAferRDfP4?>oyBLYC+=wac!#Lm7#%fQ}EX
zuy~DT+g-VkbCyyqGT2M{*P4ljPhS#B1*)Mb?-1_kA3dia`jMvLQkIYGCRh~%eUYgR
zHd-9F@0apM9;RdwdUuj?@nc{St9v>^LlpGA7wm^YY*Bqmi4m*&-3~#{zVaU6iL$AU
zp{{HB-C;&s+#|au!oz!7fOm4?RAl+DPdKy5I^R_9Xz`Cf?u8AXakhudIAj+$zIHl-
z>FG|r%io7+dw|0ftr2dD7SiinL9PKoN4MNd>xQ@KQrBpf&%e$04ctqZ_FN2YxIe}o
zTWrD_BMd}O%T;L3Z;!~rH*<+F{!C0#AL^^en?wAUYAQO5l}*t=TWK
z-??l|TwbJX4q85EBl5Mg37=Z^K};spSXq7SdePdVA5l
zipx~hhTo^Sa!2Wlu{?*(SxIR_yJDb0EW4vdo{hm@co2WaC*P)wsl3~`mQo(Jy~g+F
zfs4xNEQ{d4d{<8S9z3_aogqmIpq3a9BCCO_1Zo9?t?y^3Pq185cu*g&{g74wkrl>q
zS0Qh6lDNciv{_=fhR4#aH?2AWgg$uT9G78wf-Pk}j@}Jw7rw7S-w`z%xTsSX)ckeH
zIR{?1JLFf!&FhDKED8$rvr!NKEG$*VDVAmT)^qvA&abYE5$zg0TZZOwZ!}F+o3ojT
ztD+A+HUev*?2L`&QWIcVh`QbZU*n&pbbNS>3!b{pq8K@Am`yJ@onu#ch(!y}zhX3U
zvvtR&=oP)8!HL|hrM!@D){oo^JOUAVwiLp)6IT^idlI3MhIhDl9u5wCy;P~Tk)Qgu
zn>4w|-e!j)+jyGzCv0}e2J6%Vfp~WColNDOuXrPW=Z?pgpsYW+Tz`3xvq(&}LH^Qq
zX$1CCf8LGNsw(XX%s8W%q$R7&zTRy9HYd9Q%bX1P)fZW#S@tu^Zt9HIUeb{Eoe@lP
zVv?Nh5ch0CZ?ars?@#rPdJ6iHg_LtPS)l3lo36fol@JWErC^PDi#ZqrU=4yvH8w_4
zr7ylD(=)6!>dA_}^oMeVNf`pNy0>9V+d{_7zSiY7Hd9kuC9Sp8KqJ{+_LHq0Z&YgQ
z6mvCZiLFZ&kEieXqyQ973QQ}&Ul?;@FYaeMBb{J^NYdyJmnhh9zwN5=6fWU9MUp<%
zFvT`VHj64a3+yZt%Tq#t9=wSR+gvImyCY}yG`v`+Gn&XD1UhQDhT4n
zxQ(yhdB$i5-6s$Hxc)XV1}4VaM>~tFC&f6A`&2rcgfBXg3*2AYB`p`Do=Sw#P=W64
zL+NBja4pQn7sH-163UdKy2lUNE18EISvM3r%iv9Y}%%U=g5{=nMCdB447^=kXFtL;pd6Ow#+
zG9oeD8nXWNSy$-g)|k_)d-S0-So&1S4FJ%@8$157w#x&7CA1y4)9rf~!&=LLVmey`
zTP>p6(vDh3|A(iq42Y`jx)$kXXoR7=Q@V$4Q32@^kQk8e?hZjZ1yrO}y?X@3qdJEU%~W(FK0P!RJ*&^EK%rj|Y2Vq=em{nZolA8k^4q
zR!K<$NWs_~PS{w|v4|Z2jZ61WCi?FKJhC%Zr*SrmC3d;YX(LhmC-&N%F#@>#2&SM!nj?RLji&(s5=|0M}Z((VXk~`yo#}9Lr@ydfqO`W
zxz{sqxj_*tpLJ!c(T)&|McA)3&^&|~@OmjZ6B2SAt=_6c1n>M2Rj@1f&QQU~=$@lS
zQYz(kePX~SN~Twp{!A2d6SE6|uS-TP4_mVkL6VWNVO$JS>##!a=e6|{tHI>{&PMaQ
zMA6XMmOk%Y_KTLn)tjFTIm{~n7ip}fVpjcP>kQ9Yv7SxRksD$1k$W%3r5}rrCQpuOL+w
z9~7<1^&sGfip1aMB}Y~;u%b~p=lh&e2+%EGQtX(Yr2m2l%vc*AsPbqsf1oA2Y+bjC
zWdzD0KGA~8M>Q2ZpoRvpYARpsB30RwQ>#<+{;;PzZ(k(!LUECOGU&tH1@yr`3c7@N
z@dMKC(D*ht?)rxBpwUF<976)7V7d-*j-k3z^vq2^oby$Y>@$6jTR;CK*yfvG_|xE8W~AY(*GrF{;CSMGGrlxIGt8j^`gmI0AA-e;
z$jxjrIEYvR@P&nPmxHynYR$Zv6XK`6jFp+K{+$j3RaZVH|e;UpNXRZ(wcGZ%@fw
zE_S>$cZVgKEc{Ia-Utiu)AMV+x9VH=t>U5|8}yxOa>G#~9l-q}pu!vx%Y&78NQ=HP8@W@XbiqGwd2w(_}xwy>704NR4qFU(iT5C3G0yysE0+%NW@&koy5LN6WO0!GZ6
z7}63Rwlf`N$F{*_`W6a0hgFri4?YmP^CjTquFDR@?X}m|>0FI9vbp=CQQUJ!Kl(HA
zb5q9z)Rk$63>5>+8=py9Mi%pJoTfYs_xKhC#uCKYkt>y4krj6ahT}{Y*s$XMVgS#D
zi`33b*u=h28A$fqik$pxvDV15N$8fv((%yFyJ;He%h$-&MNv#zaJQf}xo86$BP=gX
zzK^(V+d?At*6rgp@u=LkB-(k?3Ovnz5b>#+mz~&jNIP@X`&$@szF>K)!^|jNO?byz
z)EYE4`TU|unVZw@6!qfjf%4wTA!-dkBP|2lKiOaMVg3N8RHi;k?Ml==yxH5YXlZjHdt=3*|X^D>B>kJJA1*K#=f0F!;s-0Ssa-##u
zfG$l}I|@dln~j=a-Hi0K930t89%-^T6&*XEikmbT3
zo?8coNO0%Npj^=VNV>2nKSgMgeqOjlU6b-%^I$W`@J^6Q)vYvD;VyEemIGCe1#S)O
z{uD`k!eq48@u;-aMC
z)oD!R48UVlPRqy#cz?3sCJhLJzw?m1b8tDe|A(o>Mw>gteCRdtr~LggX%+p`g7LN-p*7`|72pvOe!|k^#?YI=tF%k1G{0@P!e&5BAgBmkO*?MU7CCp7al
zC_4pF3=Ns>9p}Arf|rj18HuWy0WLZYPfbM+>qGm!KsV_LL*ji;8J_3el#@QQUr#Pa
zHjNd(V4#3y5D!0RVXtKgQEq3ChIZX#Vq<<{3FZ^+a)6>qpv2~)h?BXqT+eyYA{7_G
zgpyD0U=p`|wa4rlJH_=EE}z|E)TKY=W5~`bZ~QM>)dQ!^hFzD7ei@(!gjCw45e~ym6A#UT)llCq8z+g
zrhTMP6w-6zFP*Mc6ffs}Q#QggP#NH2`@Xl{)6rJz+TH~B9;ZQsb%^~UK!u!;)80;{
zz)tRaaW$ubNa-qG+-1(p0@+ujo1){`S-Y&muC&-c@VZzOvwa#
zl!^!WegD9+t$~J4t#V6}9MGWZ+9R&ZLfXqnIJeW|nk>j##v3I>xZt{(QkZt-Q%xp1
z6j}85X)mr(XjrgN^nX=#Wr3wiJcabgP`m-n~CR#`K5pj4G)90?>T@-M|Xl>**5FH4K<-7pGjYjS=K!eg$
zkf%D^rw+C~1JHzmZC-oS4LO|>W-#LXN;VGQwi#r}%@lKid(Pfhbp)hR3HsUa1$9
z(7^*p>84u?SrJeFBQ*#zk)-IfYX_%z*=t>C64jM0+cFWU;p
z{Zu!k{qY;`!iU;AC!%tDJOa+gS)n4{HVoO_H)HH8xO}IUa(&w+`3Oz{^p!bF@5LaN$l1G;4C;C=9wb*7|7C
zs8PJVyYU{XA@QWG;}8}bvG%ji8``4~ZX9&^%miN;tl?U6D?N3dKXB72upCG>>^Iif
zJj5t)4>T#7K$4^NLVlO%Qm+r);BUoyPK3@wXaaeQC~jI>+5$=emvMww-W1lH4x>5R(HvxqO>Js^z_iF`1r+*MiXPKJ0xPU6#+q+YCzKWf>4^`L!qg{2DP;8PP)st=33y^h}@1ugB>Gfie7pAFpoTS~T;McwZ1%->0HqMj03t
z91jMs@1@-#S7*Rr>_e{3aXUv11e48fy<$_375!>$T)i|deN(Kxa-x9{H>HR$LCfB1
zorbcf0^aS}1mtOLFiov6*m-q@c~_Y5$b`$-&0wn3=V9LgQHP
zOR3l*-p(*WO~*tAIx$vG$u4s$Z|)bvFu|6<_q^l{HkR7A1B;!5lDg~Nc0@b{p5vs#
zNz@U$9z60euNP;oEB*;A<>F|o2~>){~YdW7gpg^PMzP|f69+pEa3PtH>N2c^SrsiI4h(N(mY7!GQtSRc=^la
zSR91R@@=!<$ogy05wsy89=u=MZJhVl&b-Ir*^}$Cjh>E%&iqb=f6}Fn^a0+J_aIR2
z|8R^HeigTo8`Du`V6ZMM8(UEycc8LhkKOXd5}px%Q#Pz+-C#|v#$t2jPJCwX29&mM
zVUX{U)?QiEW;xDP^-YMxyj=bk;{53IT5HNm&3#)2G)@tm*2pE(H8wiG5o{RpNi79=j7$eq>E(4!)AjahYu-
zYkM-?{j`db-21-d{2ywaDBx+M!8i$;=#jWWKh&E5U>L?^QQnX2jkZpFHm`)lpUSUB
z+JY9?oyR*Be{Uz7xH`imRJM*P_7|zU28*5&eJX>^84w6DL@&j?loek^!c;B1ve{Q;3YEWK#5jj`0&0bIIVIZB2cJ8%9P
zo}W>Fw-e9M!%UshjQ9%R`rh`4XX>XZlei95n%YKsiPq;ZwjqIgWSv&UjWi8XPh$jr
znszp|<{Z=c{pLLbaW6q^$6=S6l*H=^hb^rkE#Jb1NsQtp1bAc_>+LKDn!v&{V;$#M
zRMELC;@jJDIebDxf{d0NKCUQkK}u7w4|cs}Q`gcVaW$q?Hjkg^)H?N+~tRDEAUa;
zg>+V$-C4U(KUBp51TO?`N@^uHV--@LyWs1*~hQt
z|0$s^*23|ybXS*qJMgz^lWa3j+bUTX?;MuCd<^UzL2}GI5UpcE5C4P`KiObXN0OoV
zHx`Z}R>vA)8t;p-#AS(R;iewWVj$mTA1?
zhZ}L46yd){oZ{HyW~)L*sZd@vHgv09sg?62zR}KLGQ77$b$^2zOS?~-nJMyCaqWXm
z0SxL{56ezPyFB&H`nX=GwRJi+V~c~mUSr}AGjJpBocFmAIA_)wfhjQXE@TROs?H*7
z&UavB5VwRQnMy?f*ON!16^>&h6=V)*lTr`7S*UriRG$Cs>&9Y!*;J+j8qG*1`=Nkp
z)P08Qv!1vt#oD-D5|j5u@73{`o9`nv&BMd`wY%%|-*d?j%txFsExSoFSVXo~D-uxo
zvr;JJ-!&cgd@YFCPHocaN{Q8ICPnRa^3Z7gFn4)DrN}4)8yi58(k%8cV{(gymI6E&
zE@@+qb*3>rTsRZ_H1?IjzNQH^t1}0q4H@tSn9@6axm?@F@}v+Ly;Y!>?<@`CrRe65
zSmTfuUBRixc*OXfl+HO`|B#;EpY8p9S2EV!0!In0K3H9VXPlCw{MS(4OQw^$vm%EQ)ynyieZ=BaAiF;w_V?M77As(_?sWZ-#yt^scXPr#9vfk
zhQNM{)V}2CvBTqkK^;R_$!r(knUvg+zk(CGm8%0+i4F_|0Aua1UI3n+-#gC39s>-F
zfw@+Qvg=Yd(=)pL+9~omnDZ6TtQYNmj`-h^o_t>s#=P}d#ARnVYDO1f&!E+ed3cAW
z&wK=%(}op@5PH(2rZ}UDv{9~>kzuD9bFjC5ba-N)dSo0Ep|K|Qg%qM1P$?|eR`(||
zQE{ZO;`v--(#T-*(_|OHW_L||r%mPR>#bY@6&73^XxtO%)ulhmhWdXFp%hlhaY0=$3Ke5Z*~f=WQg%5FswLH7UazZ#SYE0n}k_mI`t{=?l`Jxm-NN54bl_^
zuHKjOp?h}XIW_b2_CZ~!+X6waaDa@4R#vtt!EPR}xhv5fRkX!zx$jU`8J&qxcti}s
z>3y1{-pWVxO(2mm7;nyfJ9%e53SxAh_S=<2{4d_=@FL5q_{T=uvIqO4gPp%3M)A2G
z9XPUNGUVqz-Q%IQ?=0tkfGMgYNyS`p`(auO<%cqVfTBua=Y+~N#y<;irE!ekOySFs
zvK`blTbR;%71L5ph&5M_I{anK!-mmpSjvP|`$j`&$=VX0SSoCwOc=v`GmvDnK1h^6
zV)v^=bc^VgK6;DRFn==gfT%;_9xiN>J(E-@iEl8xrS5Tb_WW8M?O5^LriObFyh<9f
zP;`#>qkx;?d1YlmCA-;XKlbMwlVv9C`u}AC(v>!&YH-=h!PbH{&BOp{t<99U{wB!)Fax(87uRNl4iufGADuOpdg
zv;@210Jf?UPfX)a_qGjoy6qvjB0RHIOjl*xAqxYw%~S5P&M#ULp5jjOvHpmPHCi4o
z;uk;A#630+W%z#1*|weguMm;mO_rKto?X&if5H8)Hzd+0+XQv0Z#-E>b096eyogtJ17C
zjJJlNB|&0~8)(zv@Qg1e_?0i||LK|DKK7x@
zd5L#*jp&nMqH4mUyd-<^J|BWF1d_;ttf>RjeA+2TUAL8(CbHuRkB;L7U-3Mvln34T
zE7M#o=(MjZ#?EWeUib;Pbv{5>vsnMkou0==v0|}CNFFfFz3sEJ^l018FBU$Z4PcEn
zPPP$qtZAa4Cm=1}gOR5QmDepV1CsAM_Lc_H)QHg^4v3s$);jU%yT`;oh^6Rx1?^bH
zcn@!J{<=9uM4HtXDUG25xHsK$pCB@4%dQKBFD&QY=on|b(OxLvwfl^jg`HeUv@3Z}
zf5J>AJSb^-NtghQVc23J=XzyjY}CoZkrG$)QdmICpTV4m0U~q%;p$&4A78*C4b9zl7%T#oOjvN#ZfX>;y)M}_GKTOd-&>{
ziqgJ48{@)8qNa0B9F(WH+5;A!#hVs?NjVNV&4)6fQJKwf?hS=BHAvpai7f5d`
zGs}O_gKEshQVM)RPHncm)694i+CVOo`BIzT0d%Xfm@|H(F*huBXnRbt0+t$psvPej
zUGIr<{6}yj&Y2bk*GR!=>APbj_p
z+-RVpe8ZFrH-5dY#r_RUuo8Fw(N!_iTxup@Nn>g&P|IiijzY)u29zSCGA@+8l~*W2
zH~$i`_jNChzfPy;C4-t$|Hclbat!&f5Ogxkl60JeCm=Ho{dOes3AvOQO*c20AnpE
zx$CtNVFa`wDOz|th0RFyC+bZ3U)70(4k0G2V1HFXXf8#ZZ;9oBlW5dA3-k14qY5MJOb7;;mpmQ)E-Oz7h^<{rqJbmk~r{`GD=D&)f4QvHt+&9?IH
zCT9I-DPr$OY5WSnCvQ}r{H^S8Geg@G!5!H^?01W%s4EJ>HI46Q=fjn6w$wpX4x~p5
zc%!I(IBs=~c*z`t`vZl6cqKn=yxQm)Gm@By$TqZ+rcDQi{|g9@Hu4WyI2j
zX}M+|f`z>X_vpc`rV8%Hl#b^ur}6h^3vs@aXmQ&TDAL)w9xu5Mwo!OH*DBvc=sH6$
zt(e^q-My-gk^I-hmfWLH=2gLU^uH`hbURM%RAKx{otHx0}>mfIN$j!rG02pB}6Z&Y5hH8
zA9^mFYq1Ar6UHiE!$i>QnjIF!O~uVJh3_0
z+#u}BA{40kz7GwDeg|HKx-Z7sFKy|{ALShrGkC9de0L~6Mck=}O1RRVy$6>$=hFN@&}b%Ou3q;D)8
zK_MPVWK)53gzozK=w*YM^Jt2oYm?jEk*j!
zHdw>61f{X-VO*&sjV+ZG?}Dh=`he(#3CUcI@$2`fJ10-sDUKA!FSeIE70z53DN@6w
zC)>3HgG>@fG^XH=mIM7dt5Rz*r{U%#Ah$#vZqNe8kkSR#{4)mko8pWOF8M?9+7h-}
z493oO^r8SGMX5qv6mT~MOLmV4ciOGJsP%UFQ6}00Y3Eq^!jAq{=K*4{U}qaw1HJkr
z`rv>(`u#G3aBs9nGUNQjAdTv>6-QXzk4gpiA??ZtBX%gm&UVCAf@j&B=&Ef1vLdnl
z&LjG;mHXbebsx&9lEyjOKNjrY^G9um1?Q%QkiHeBa5;D6VzQ$$F$31NT3G9n=J$^N
zl|g5jYz&7@d3Cl}k#!tMVA%g^3*9Sty#KpJ%5_=cMe>deD3rZuxFbF~P3`(O+
z>K2TaVa+H_d5jd*ojc$XBfEisZQIu?&A1IeCUfX}FJTWO*V(Du3s)q(Ms3?^`&v#s
z=Jlgyb%caS=j7JNpB&tEBck1#IpQ)Z`Eg(ZR1wQJ_#Cq<01
zA~6LspGjZ}lWIh8;%1xTuOuONshNVB9@0NElddPd66bpHCrwuK-;v8U2O_+FHky;3
z-m$B>sVfb9Hl@2GSnslj5N-^`+-HTmfY=8ifqct>601-f>WC9j6UODuF)K=%6kZ}>
zG9O%@ys=w?g-_g}xBvrS#?HHtq3a_+a)-DN%@axdSF&U$@_kx^)XPbMOKvL>Hux%~
zdmNH)fwwC97(j~TLHZHi$+i3xyZ8uGfP>8sQua7Q{}M?Y?LXzLO`=)kOL)(ys#;5S
zG8Fe*`Fv$WZkThKh)MITsMUqh`(pDJI0LU*Etr`w6Z;4SL
ztazQv%B_;mi;37BGCUS?0Qb1SQq(8nijz(l2w@mH-Toa_!dC&}4;JEDhDP?BCJv-*
z2X`j@7&M635E@(RB<3_F@q1;1yW-EIc1TQTA48z5xKEA*!&l%gqUA)=)wSGsAci>u
z=LIB5A0F^_6l}tN3k@Zy&YsGFDsHO@4jr94^;-3skm6m0LEB2_J;?^TK#Qfsj8AiBg9#mhL47zzpnPkUaU$B^BF*
z8)iXMY-cYLZn{!O7am`=90ajKp>8>C{V&};-9Y~?$UgF{B!_V7wce?cyX^gK7nCz|
zK^23<;j$@81ZxKcv7dQduzXKMf`ZEf?dgpt@vwgWNS!=SR$4p`%|7En!W(53=NtvC
z%S*HAz54t(jd7{(SALL?7IXmJwOdT!f7#P6uTU-ePB1@Fx!zKF!lA?CHH>{t>X!&j
zM(F1kC&HPjE@_{kbUR-gFuaYwBrkby_e#i~
zGG-6-vgkEwPZXPK90=rPsWP@ZC=?{_N{f77g;(B1G}|v1$=sPs&IlkTseMe6ISE5O
z=+0p)ejt{d*`uA}U3+6MoQNWK0IZP1m-|B^0SMHF*7sj88iVMv8vs%$yczY66QN$*V#_ADO0h9}#MJBsj><2jqk
z`>#;lCNCsbX-{kkC(!rxKW(^?inGdsD4B1zo5S3JSM3r=;p^N#L`YX-8*Thx<|xn;
z1QjbfgfeW6SYYuVKJop&J?MHGxEV`Wj!g|HmApUPhi(H@Z0U}FZMMBlm6NNJc-;Po
z%6xfIR<2HA7e=3zp)gXR0HYFhtcXn&Q>cOrN|b!MV*y#cGlp+REr!2Vaw5vLvAnEf
z4Ti*b0a@jKFYVu2_()z-%WE<`B>YY{I@y#O|8bd_NAb|^Uy-Bot^#@jc^;q1LZ4{5
zzTGPCinh+r%b*z9Z7n@*cD-=z_>izLz(%qlIckeRgBYjP?63M_{8*|g69?$M7;w2<
z{6Vm}BKr!fgY29m3F&u5&i#l5Lm$trYd0jBsLN_@`o{}Hy6ve9qJ#(X)tSKf6I0wV
zN>W{kyW{CF1O>$5-1LQ`;nw;mPQmtZ|B~=dKXIC4lT+wyv@O-o^sD}3K9Q7awzzH5
zP3e-nelV3nRP446y;#0_yTD$Kd5;{KW3~)6;dh~>uUNd;$;5^XiCJUG>H~^hXhYaw
z>Xz|O$>z65`_wq%oB;&Kaqv}U4$&tkw|Jf#?+(BxmmhO)_L+XjSalQnSiezcx4_c`
zyO~ssUQ9hU5Z)&Xne_d}imUe)`_e|o9sjIfkf&~h{Oh*;SzhuKOynz((!BeG{95M(
z1Mr$-e|}3c-3h;$j6g>K&k_l^WCK4kPfT|d)$wTdqTjV~nRffMzg4An2%
zib7NL%FKON?vEOePn*ze3Ld3&EGchBXJvU$G+~#E;H_#woL{aI3H5xIq8rdDc;kx~
zWqBo^^`K?y#1!e}3j&Mn=d^n30VwDcpj*Y@oWf-ZqVZmgrn#G=}W|h<;>p!KDif>Q*
zYCHwsO?x`R3&V(HQR-VDV2McE(5g*ao@S)0R*_Oi?;26`MJQFDesX<-X6Cp?#Jdxy
zY)b1|*G4zt!fC7GW#ih2bbvfX@C(L^MmW*MSu2vAjY$j+RvA70lk$f2r1-f<2}KKO
z`*8oaum`77&z}_ckAd|AFC(MdC$A+J_U^zj0sUG+e{(C)Oyg$zr%vKT^lGI~8
zMFlg28H%;`uhg6-WT>ZUq93Km?aWPrQz(tY;io64$5C-|WJ><;(EVK>!rl=nbMVJL
zyEb<;GO6~Xa+jEis3RINuxJoc9w1Ms&ys=J6Ppb;X;#+Y`f9G}xwKPyyGw`pk|vfn
z1Lbq`wq4RFIMxa0gapcRlAE*?thuWC4p{5g*N_6KwZH5(C3h^q2o^-$fK*T2r;BrByxIRLbx6X=KU>lEU~b@i~||fZH;0;En<-sX>AD`j%2zB}+!>q7Ag}4U5`8rWWJIfw~6LM
z)_FviVlmlYmB^0*-ECe!3CeNJ>I-OAP=eIu5Awx2arH1ol{65vGflnRbPP@E3QZbp
z3r2D7?5|NI&bT$gYpuQ|Jl&gw7VBL~&lhoQx
zGsH}G;h+rxv%i~NSE(Ama&QZ5RQ-B1x9_?5w;2avogRGpiJ6Stfmj^8KhJ)x6H#C7
zZJSQ*mLezz82eS!^0|~qJUFklq=O|{%!Rn!DPQKI!3!FD@ABcQ2_PmsMf0^)T)#Vwy?qxmT
zvflqD*wuOaKvZ;Lmee60K6V(V&!@2^Xi-nYndp0(C+k`k>%4HOdx-f^_4gF=TmD`t
zucSA*YUB>UPawY9n_-#B)xHHxyNAHLC0c&>RhCH+dfWAUyZQixcu=2uuD5^1Yb5^X
zEpuY_p5So*;EhmQ6i^+0`=8j{wNCwymf*38yr8f!{q(qM!<2A3NSand*?=}
zF5F8|M!Hq`pm0y=nH@B2wW8;WknE@1I{RyQDSl2=3dxJ$X7l<$u@vT-M$&@4TGiPD
z(bMA5&oxXwq_LxKMuYsRTM}Ct+E*)bB=z{cLfY(DUu}3(u{{lk_~WrgEEL==a}OCl
zK5bB0IH*4v>*m_CA+G(&59-~xed2Fl?|p9*)2k>S|L1!l6~}&M4V$g~tLm;I(k>zu
ztOMT4LLiFn?I6DRb(g9{*LHuK!ti6ykK>7v%iOVQZ{5U(i5?sY%&AQOEZ=<7m+VwZ
zJ!B{Cj6kp_Q_rEcL|Bu`=h9k$9t5xEWK8GgV2_~ti~lO
z80cG|17?fkw%Kf`GDFrL>Pt+mc}eF!DGbFifPc!89d)dIq2QRtSQGqDUAXgU_zX%E
zCnSEFpm+BpreKL*k>dVvguUN!zkJ9{RBv6+Q~6atteHCY@d?k0Rp{dW#Y)pJK}x71
zX!fA1-6X&E%!gQEWAC{^NC#5Y`N%2>P+Z7F7VGa!<9+g@5f9PaCE7{@^=sLuh+r2h
z+Wkp2srg`~gS>q80`O7%{@rruR$nuXbrxv!II^`7FAJVF_vF~)x}voF`CRrYPc!=d
zO6PoByO>tfMt_Z
zq7?Y=jaCFvip*RUj`T<%#j%qAb=h2
z73bTyO}*m25pPQkLpcBdc?v;h2Cj+7Z5J$7X^7FOfi-j@erJLKSHnAja(dUkZ03jm
zgGkCNP6mzlmPcz8*6-=$?rB`0&TsD+D?-%L?Vyi_s`CI|@<%)+M1pgeBL1e__-_mT
z-P|Bnyo-b5+YF(#oAP*s1lGUC+}2#Yh2?i6Ij8E&M-ee-InurJ^KsXTxRtT|;>$b*
z^H+`gbqcc)Y&?YPOM=VjuH
z6E`;XlagdM^-*9VTODuq3%`{faZ;Ab2x?}!eB2K?Hd6#xI)>7}xx9|N^0yAs~00)Wr3d%tF>BVdX+OaOXa}yGd
zOOX@>1n`YNA8u;+>zrY;PiPR>FPiKTGpbb%Y0UIa|DRU=j-fPc+B=X*($8t#&`YfO
zAQwg@ch(EewhX7@-`k8t!G@bdo8Pe;g=lPy>Q4ZO4x7KqVZ^mfS@Pj#n4YMTQ~wdU
zS#(93gjG6$iuAaH2ee+FwDFxnfBN_LDR^Dd)Gr<#vVZlzXw@(hIg2m)O#tqm`yaL&
z*CJIRco*9Do+E&YCxYD8D5sAlIE@1uR>I?Kl!HEQzHb-bV`r2p-%P7bCtvBLeHJnv
z(Ay-P3EVHvW3S;xe50kVx%6`TvU%Nu?qD;|kRIRg3HF5jZiOm54&8>fYW3Gt7+i4&
z=tV={yP#6;DK7|U$M;@1(X16d=jAVgVf*68a1wm)6UjopwN%Tk)wNMLX^ZN4^>C7I
z+bX*3B)1a*WjX*OeC2bqE0~bbO{?QA(}C-x#6QL|ur5_}T;77>9a2iqU@*>XNfIel
zaLWtr+*qW1i4;R*b`tZviQ~6Dec@YSzB&%EIAtB^dlqT7jwbo;FL6_oS}GA;FT$`{
z-;l(G_45)@)K-|P?LUr>_j{*EsovRU8Se*FI>t$*KIZS~^gHB`I=y?oF7n8e;P*tw
z&l!k2W{;;xsM=c0H!G}atoU?!9`vD_Im=SLU??81H`$s3^`DUA9jjoWSW){RMvpe4
zfBUUR5H#zc)BT?2^6wm^p6IAR`+YW!RJ`26Js$a@-WQ@D^=jC$3;iZTJ?=al0djXHdEyMl?kz-CD=Gzi*unXZJtrTW=o*>945
z#-DkhbX0?God*!7i>i}TT-_!?xSks}XUxe}V4?f^i_{HxaJIFVM`oHM*Xws=xmBHH
ze?NNvJ*mXg8@XlIt7x5m8bd*lCpN7nD#@PW@z!oQr!wL3;Z1!f)<-yEHtSmCG2i7K
z_kK)4cWjoYKK^{afs?t(?b0@>)5tbnaXKj6RsZNBK*
zR1;R$L|!J{Uf&hA4A$FG%P$GRC?Su=u?=r3CG!a
zDal&aWT%7?Ytujr{^5ZeqPw!1DxB7
zeHrxLvjJewB<(uotN>xl{NRjQlzwwJ-=A$5Id9KaT&!&CQJ1ZHY?PUQBi!4CTFga@
z+qP2itxch3QI%BQEi-e?p9tIZh1LeCDo^3rAX#$v$XaF{`*D*s?NdOwdn2I3sM~!%
zD^O6zA5qP9Ba6a^=*U-e_};@g$Jy0nb19JQ;*t*twJ9>OQZHp-ny|ko=X%CnB6{5P
zlfYjXQbH!0AR&|Ae|<@MbqZIhG!1T>+=iS?iFOQDNu#UWE7F7q3<_wwP5sz!sd5-8yD6=+bM4PVwt|n6%bC&;TY0X)Nxy_4Y77fL4DkfgshRA
zh|fVX6h~kDgbVqQK`Cd{%1M~K23;)60IeNm=XR)v)gX7ZA!g`aVfaH{6xDTm!H5@=
zhDE$V{Qt55(Kr4m{>Z6v6SWcrBlR|zN)D~@9!v|Oc-KdFe=&9}*NS~{8p%|6SimT5
zCo<2Gv1^^*1k^Y0)DSC;UmkJqw(1F3G^061~W>CvxT}^!v0J}bBAwYr+ynt|p#yck`Ai_x8AV8G&B6TuUx2Q8iVZoO!#YXjB
z5iH-Nv_qc4Bf?$i@YZl6B~R*YL6wHZ8_@nEIZM2oX}{y`aGxVwewQK@oko5EE;78f
zgM=7d%w)4La@}1{|4jwO5tq_~SL>W<2HH+21b;3@ctk3El0`%~o-ao;OtofyGOFNI
z{{;(0{RIwFi@Fx1e?pTrDd&nZyC+(;gyNY16+h%`pQzElFQertI-2SUG;x4iy`HEO
zy=L6>ssd*e@Gct+G0fi3WV$>WtuS8Hr}6KiUV~=W3T?X&ePv0mRkI@MDZ#RNAN1&o
zR#;Ro{vx~=>E5Wr^ARs0?O5quF`nP&wY?fPa=QiD&L0%CyU$G(mfarV2bS!r+Vx`e
zMJ9|J&rP+RI{*iBtEa;5Zp?MIO_`18F-ByUT~x3cSz6%yUg?^Zig%cyl5cx}8`=8N3a4!k6zCX}DfLaW}3L_8Xs`sJ*)~czabp
z)|TR*m_Z*qD9XCiCjXJ1eWv7?_M}lDnHr*A6}K9?H>9~$4Oq=n&>qERKmAlsL+|BS
zBo#r5Gw8MFAdeWgsGL0D7@ICJzUvLMYYeD8NC1ftUhhpsz>(+^xCbQ%)wBt_ZpE*t
z<|0vG*w{_dMF$BufcN^M#obcmKNpD$KU3Xy26IVrokBg(w
zU%&Dt1lyXy03`8V9Jgr463y<~<3;fXveOM&7K{zl?dDhbdp(!)_13{a;)ziRJ;ueLbrK#J^%xDiQJ
z9~yC6JZ~P=tDo*P{vSBT3*cfIrQQlue}`%0vVXh4g52uiD?^TonI3~}h;efTV)iAD
zcknEm3~ZP}_VDa_Zx`ukiNdozl32pV5SsJAKN4R|&~!-<{D*LuXuEdVM4_!<%S2_S
z_+~!#lOf;o(?whh;<%?QBk#>E%&VOL6zNr%;bluuDA&T*sUjX^-d~DcX!xo3o_iYj
zhvY23c@M_($;d^r+Dhay`*euFoI0E31{`jfUQ18Eyl(dPgEQ3Egj|nRHL1
zZ9WO|(x41ID0T9FD}b1u%;3uv)2)OE2mo`R*)0@n>|3ig@V0-_i|neAE&=SM5!29@by#PI=*>lF^WFXcbH(y#5N9^+pfvTY%Knx$^GfnCn#lUFKddb_cF697*TeXMlUU&Od&
zI-KDt_IFotX`@KCFlc7?824!AGt0Bz%mkyt7B-kMd&85o?P!t90Pr
z5#?fhlu9^xvu4tMql$#=K)1vQ1x<%wUuHl5OE%?kCbiCMH{#)>aETDl7yRnEGB?Yd
zUw>$A>>?XXt2hpbXjbkf|8!`V3DtFl}VIhP`xz*tg)?yQ&T*UQ7ePR}E
zahx^NNmJQg`sb){zVQ+`{n8ovh<5mr_HSb5x_d&9W}O;kE!}yX%80PIIC65(F3M~<
z?g}-T>5`4ZMjL%p-kL5na$>BH>Gd#uTvLH>(TZ#y>4OE)hDuc>XD(!n?G*LX&Zqhh
z_24i125K1fc5;#me~TMICibnBz2ys(LRIRoip12W7wGb5X_%cr{$j%St=w!0t3!C}
z#>|I<%QX+XpBB*KHyFvg8?8;cRZtj&6CKMRmYO}kvO$?sqijE(-?f+H^GtJFSlUya
zyUiyp3t|?*i4uO00M(`&jhV#K%dYbk&NSuoV!U^p-Wf9gup&fp={3uZP~~W1^GSz+
zFe@*X)v}OJqI`K^%68vr7KJxjN#>WB{V0srdu`*NWH$ztNsoPUx)~8a?;?q~6;~0f
zn%@@suv8wSVDg3fDF=gRY(?bGwL4~uEa>5B-KJCw;v{|hVcy9QCyox=7v10wB&BP^YqsoV
zSbbFtgCW%7Or_*O6S!lSi=yUp{yf32)2_>-Y+JcK5-Trt?=Q+$jH8D&4~jeOVvU{U
z=@TqHho3$s4bgD~$9J&E=*|cD~#&$w4>EeMJQ?=3~WmRld0mkUW*N+wK`7wK(7N
z=HI4s8#8q%e^>yWXbH}CP<%}p_&=uJGA!!u`yN)LQ@XpmyQD!tK#&p%>FyX(>F$mp
z1r%wdySt?ugdqop7;4}@+~4PSUC-N@cl)f*-utY*_Bu^wW|p~X3fl61MDcy=#n8pI
zbZNztT$e}nnOaNycuI9%$3=nQ#o_hpr04bbcJ|m#n0N58e~2d;Thvhp^RHbOfe&3h
z;>PlZLlKNPvL&06HXV1oD{1{_SQdWUD>o#q@y^BIOYy`)(@ghk
zRuQ@T>ghE$E+kh*K0F-px*e?HX$Ndta%*B|$;rX>IkaC1U)|abF
zsrAo|UqZ1fE7g~pH)kygDZ`HMaz+kXP0?K9ZOQ}qDaTH1t@gcg6L{@a{Yi~N|0}-p@&6>z!HM|~iq_~JhdIT$Dhs%+s
zhDNpMqZSaj+3f8Rzs}?MDjF376oj(@1jqp)kAPk3#3y
zdy<<;h{J`-m14{=50ieVXFtkoa*-dC}r#q
z0Ti`-Hg&idKXTu}IuXKJENLe&*v8%*8C~Vv`EVO%lGf@CzG177%)vAyg0!~7mFra2
z8$T(t|DLD?hke>mhBqK%uY(dcs{bbT?5nhP_#PC8J^1g@dcFEk-6~gITVQ(&%ynd3
zKKe$Av;)4WsX=kV-9(fM<8S{PhC|<{og2AAs8rB_c%oxYjXXfHCOpKyTLvZ0czp>Ul^kJqzAW)fc0q_&
zuB4mU$Zke04QyndemqQQzMF4i2k8-#$bBGG+;NIO_9pMp9CF>AI^0vu!6fHTnk!JTQOTD>I0o)+OcmJh0dLS1
zO~IK1;6!9q3_hM@;a#7x@rR-LVxPU6>X}nUKQs^+y6$SZIU`CdBO@w(!@Rq{>o_c)
z6uPeUTvwEqJxEf+9Y4P5(0#9~gzNE4JONV=THhqn?jV>LtQC?ng5P(P{{IJ%9(r-W
zI0=xPWwbZ%9sE3f-7l$pe_`_)Kt&bJ^i73I)HVvZ5Ot%W%Hl3_1uZttsuLj;j%q_0
z!7}<*_3XS+n?YOs2iwJp`t;UE*ya6vJOGzdKBl+GSQB15y#{f&W0Le5(zthDr!4}*
zAl!u;{VtZtLk2}+7-gN(#nS3ML^6s*e_01`BXCO$d;?auF=?DfSQZh5ih`>7=v#@3
zk&1Eqp+&XFi=cH0EK&frJVu^ll6;ItG!A)b+wDd8%Cp`4lL?vKb89(F)$Z4(!0hRI
ztTU9?r)Y|)2IXuSfE$M!Vj8d}xg5~ln%&(}4#pX}$uX<>AooxI>e;*CE)6$kaaxk<
zK;^ikCnI%*vFwCrdkTLg*tIneiHPQ~L8w>2$Jw@(UNkd4hdpEFO}yX1Q`R2ut&l*_
z2F3*csoztfUF-9+S&Np~u#3?jVYi$Y>C!jc|BB|Xp^Bf2^)^e*I0qA^RKYNx1_s)w
zEp^qa298NpUgJTH^RJepUnwj{C+I@Pcl(i9dePzojrD+G#5e_~a`Sr{MYbu)_{m
zk_UQ*C+zfV{k(X3JFZEP20lsI15hU5BZbNm-g-{Mfe$7rxtc5!mD`#?<6sif0}cIf
z6wBIunexg%Y=$lfqM-{{k^aX#KEw+I`CLggZkb2QhF6lW!ah64Vd{pw?YF8Ko;~}G
zZ@m+1rrzyw*7Z>vwoiM5b`|dwjT!+ig4Q@{y|SzkJG-a)TH~BaN=IE;V4~iQT7E0h
zkNHwTt+)20tr~LTQP)viHnJ2-JZD{$<_D`WQYneBWdqn!f(o4cNrZv
zBbNjN!0U4vc|5hK6<(|lik*eZynnBB@vB_Xeoj}`}ke8edD-Nu*A23+a^Wc+j*cn
z?)SGrY-WBiuswUH;T=JKMod(bG^^(S)DXj7c1z;UjN|mKS-$aT`GR@
zW>6FRjs8eZM>8o)Me^|a`JJDj^XPwS(ChAD#m!L(oTLc9BD4;0TFd_w_Gpg3m8xT9
zviNKrnK(VX`f$+|uYl^J|BHHJ9>~`CJ?Ult#n=@_1>Dnv??Azlin^v!+JDVdIQ
zTk9C!o8fUkBfFOoaO3h{#mgdIhz^r@QPM@i$G_g6@cS?AyGA`S878U2*Q#*MJxT38
z#2-<8s-7=8m)sx%JU^@VTn6E$IlqqpXg(5I4^PO!tJ#2u{UEv6@EblXI4-*b5)y`4
zh1Dd!-Ftjg%xVEi3RJqT#*c_F#}}1s8LjzHcyrj$p8y|!9meloFZ#b2phDCC`
z;5EKm)HA3=upNctW4~X_A@_S1f){Kqn8$wQz9UFgk+7GAiF6Px{$Iy@3R+)ztIo%*
zv|%z-tngPbQr?K3?HN54NvP96GLUJL*C;rL!L7k2OMP0Nabtg<^BYaoxc
zcD1mbG3oKtn9JV4f%g3z!{dPs$wvw=iONdCLnc2lkd%pfn+F9JN|BV*vWL=rw@%6R
zouaAURtgKgh}QP^5JmG0&hRP-O8)uxev9l-C<=f
z;*R69vi4_U+6r*gXlLXZRJ3H1UG+X)DBqF%X0W4LU0#Er2kzzHKh^uplG%AQRUUx9@SEVS718~sJ{;Jg4F$riZnt_~3ra2+6R
zTZr#`RZLqWM`dLzs`2#BGLg*UU(Fy&{iJQiaQ2U-YC+lZtEi({B0Kw94XDob8z8On
z>I8xHB-*~=pTYubI1faRo{^4JRDbT|u{;xYBfsa_mv=dmpmc^7eX^o6_dr<>nc$nO
zcXK3d8r{(@f99T%?94KA39|`WuM@3~K3d{pt4Ejx+;o)3C?7M>Rz+Qo=a-l6XB};iZVevTVwf~YSpO>y4scD{TTYhZ4aseMYYQhsc%zjJ@0H>SP`_l
zpav`v3MV{?g-n^NX)o=QJE}#05&UmJnwc07aS~hG2KW~k!z6rzkg$6Nc}z@dwBNG$
zPa1cdC$c}O3a`#7>p4O(dYj0EtDP5KZ|#QKJ840$HGZn`Uf;Kmxz(65!c^+#UA
zAPV0>>*-&x;}9>Q2(}WzLD<29%6Tx~>s9Z4j}`U0yBf6ZFNE`-7y&>bC=L(YTzIlC6{qYt}q}TCa&Q+i`OPG^Bad
z4+df0+=iAh7j)XT;YGu$L
zxo(#vK>sDPX6rF-Z{$h#YbFopT}Q__RJk^swz3AZ&Bv|+SCASwSKPH9vd|3Jb&Oqj
zIJ+dw5n9V&&_^_r$mkD}f9T>Z2|NkSnbS*3oC8*3Ctu{Ays(rNwG?MzYICuO#6)gr
z%Wq)O&1lm|gx9rb|Iux}?bzyM#3{fNdp6(ay5v_5!LmQ*WPipa78_XJl<5=g3nUPc
zzacIzp-BYKd~681kk{Ch8<8!wOs(a8B1BrfITj6d!rD|Ao}Zh{(F6@@QLEDhtfyWn5osAv&e{_Hm^5cJJ~jK19&Em8u4yDO(}%A_81Mwko{*F
zZ9hj2R*NT|^#AOZ0kwV9mGdT}Gatv@U5^CakwC(w5yJL>-!KU`6$&BJGEUBTOnpC{M`tlAU%PKaF%tbyW+70y>DSkA%Xnh#y(%UCFd*u
z;1kr9vHf<9Ke6e?xO??*pJBAo()@2FtD_)VO!(8{(o?nK0}k0&la1&cIh^Q9mnp5;
zC6VQqiIQ#~U$WO*gCs=hR7>OUVlamTM|9}_boV+jylaRzw@bBxNg4e^Gi~08%lJ?5
zReS(h>XsxCAYCF-$S9XTF*B^5$9p+Y#?dl>65E*g6rv$)_2<~u9Ef5S(%GB5VUTy(TMM3&WR1{
zdNmyZ#SoS37V5!?f2Bk50epoM3O2E%$UPlQbaUMD$KUviZ`9F$sD0kri%Mr`_{&hb
zjd^@gZu?-|ZM+TUJGX^m-nwG|3WEGg<{>ZEjATqJv?>)**U3k8T^x9&}P$ELq;?;aTrjXPAR
zIIR_6{%mE`Fw?B^j^7tcqf4-CWfp#7_4}>}i7tm7Qf3maqFmvuu;ITd=VZnqmxhNl
zz2?>9{MjJ@FX|Z&rElQk&IXky8f+8qb3jAp16K{oMB?l*;B=*}PU#JdC3;NXx7`-A
zy6_9e!#~$47m47n>J}U#2%JPgA?-K75^)2gWk-aGMH6}O1V5wF_GZ2wg_;bFQsPjP?e
z6<%<0Jo&U>`d{%%RVmEWW*0PH>U_Dug3Iqy$pO!>8?XM##ocB>F)#Jr+B8?iBKE{>
zw@gH+H{KCU;BIil=W~JgRFT6V#=M{1lDmWmi^u4C`Y;x>Y{!lMV8;QUE{O)`J8Sf_%*Fcv^0T&^M&B`?ko8k
ziICt*X`yqXBmcCB(<)!&Nu8s?Ub*Qiq}Fj@Nr<_zl44bDXMZbs5pEw+diuTS{iG7Z
zG^_GdfHVPTE#Yj)Swv_~G<}Nn!qA;6yT9ayNzvwduBFtfX;-u<72lJ75Vayzf}TRq
z5Pi_NsT(~XY61cvFG6YBlk&?>uGvUL8Y{8$)WR`kYfu@X>*jBUs7hz0A>f2#{UogJ-vfYWuMuFu7
zx&kNfJm0pv7dIa46JEKPdrVR?PW~io+ZyKSQ*0jTOtLsA+_X*X*CO86VUrxnG$>-
zE{Kc&jJVZmo(--z>P>OZUUoduxVs8dPJYok&~g7qqril@l=3sM+S55>!Dk(|W8Qq>
ztHby=yB8O<(eVxAE(mzz@4T}_FtWZk0yv4Q+LjzvZea0x6kV~F;KZP>xs0t;AO2rl
zy=x5DSA+P@?$m-0pMcpjx_O9DKTST#n6sSkxl?tq9|Av@W<4^T-37+_GB8Fjr07ga
z#O$6L>|MrS20ik~olAex*~D(PfAmm*%M##BkB#I)cSmeQS?V3>0H=f^3rArGexAgC
zHInRi=$soYBlw`5m95Lc;dv&qdrXcJ5u9!Vm;Ip%Rc~Ui0}=Z|Ik)noo~V>IS}6B;
zH3vvKa7R%cVoDW2fr4S)yF?t6dxi52x<}azehZl7+FBrM6qxZRhNP_UD;Z~*tm*3H
zSYZDDX#r%T3}6ZyMpx}HHa2(127@f0&w)j~w_j09yzk76KGM@u-A`nns;a`^+XK61
zkQmKNddb0fr?q?z;@+@$Vpfmlno
z-sK$OMpf9?yXs=>+X%yAbYFSJ*x%?3#X>tsD7Fab8v1fz-muziNuuONAJm=HGUp04
z`XX^@02+p_cKy_f07AG5yzoW7Fh4to&q3jq*`sgi!gh=
z77_ImL<;Fhn=BKBDXZEE$p|8X?S|)$I2f!U6u#1s&bHDqvDc5}e{D`qBT$Xos(9te
zD0cJQKb?Siw;Tkd7<%eVpMy^n@Rbn;hL;%O?k*8-U4OuLsSWA{VRa^*zc`K+P~n2G
zWy