diff --git a/docs/articles/areal-weighted-interpolation.html b/docs/articles/areal-weighted-interpolation.html index 025edfd..b007076 100644 --- a/docs/articles/areal-weighted-interpolation.html +++ b/docs/articles/areal-weighted-interpolation.html @@ -30,7 +30,7 @@
@@ -89,7 +89,7 @@vignettes/areal-weighted-interpolation.Rmd
areal-weighted-interpolation.Rmd
The example above is a spatially extensive interpolation because it involves count data. In areal
, these estimates are obtained using the aw_interpolate()
function:
aw_interpolate(wards, tid = WARD, source = race, sid = GEOID,
weight = "sum", output = "tibble", extensive = "TOTAL_E")
-#> # A tibble: 28 x 2
-#> WARD TOTAL_E
-#> <int> <dbl>
-#> 1 1 7992.
-#> 2 2 12145.
-#> 3 3 7344.
-#> 4 4 8458.
-#> 5 5 8783.
-#> 6 6 14050.
-#> 7 7 15840.
-#> 8 8 12188.
-#> 9 9 14217.
-#> 10 10 11239.
+#> # A tibble: 28 x 4
+#> OBJECTID WARD AREA TOTAL_E
+#> <dbl> <int> <dbl> <dbl>
+#> 1 1 1 46138761. 7992.
+#> 2 2 2 267817711. 12145.
+#> 3 3 3 66291644. 7344.
+#> 4 4 4 53210707. 8458.
+#> 5 5 5 60462396. 8783.
+#> 6 6 6 64337271. 14050.
+#> 7 7 7 101268146. 15840.
+#> 8 8 8 45966410. 12188.
+#> 9 9 9 73993891. 14217.
+#> 10 10 10 62915358. 11239.
#> # … with 18 more rows
For spatially extensive interpolations, a list of variable names should be supplied for the argument extensive
. This can be a single variable name, such as in the example above, or a vector of variable names:
aw_interpolate(wards, tid = WARD, source = race, sid = GEOID,
weight = "sum", output = "tibble",
extensive = c("TOTAL_E", "WHITE_E", "BLACK_E"))
-#> # A tibble: 28 x 4
-#> WARD TOTAL_E WHITE_E BLACK_E
-#> <int> <dbl> <dbl> <dbl>
-#> 1 1 7992. 153. 7779.
-#> 2 2 12145. 1323. 10639.
-#> 3 3 7344. 591. 6635.
-#> 4 4 8458. 160. 8203.
-#> 5 5 8783. 1526. 7056.
-#> 6 6 14050. 5840. 7439.
-#> 7 7 15840. 8220. 6629.
-#> 8 8 12188. 7604. 3796.
-#> 9 9 14217. 6838. 6413.
-#> 10 10 11239. 8703. 1667.
+#> # A tibble: 28 x 6
+#> OBJECTID WARD AREA TOTAL_E WHITE_E BLACK_E
+#> <dbl> <int> <dbl> <dbl> <dbl> <dbl>
+#> 1 1 1 46138761. 7992. 153. 7779.
+#> 2 2 2 267817711. 12145. 1323. 10639.
+#> 3 3 3 66291644. 7344. 591. 6635.
+#> 4 4 4 53210707. 8458. 160. 8203.
+#> 5 5 5 60462396. 8783. 1526. 7056.
+#> 6 6 6 64337271. 14050. 5840. 7439.
+#> 7 7 7 101268146. 15840. 8220. 6629.
+#> 8 8 8 45966410. 12188. 7604. 3796.
+#> 9 9 9 73993891. 14217. 6838. 6413.
+#> 10 10 10 62915358. 11239. 8703. 1667.
#> # … with 18 more rows
This ability is a key feature of areal
- iteration is built into the package by default, eliminating the need for repeated table joins after interpolations are calculated.
Spatially intensive operations are used when the data to be interpolated are a ratio. An example of these data can be found in ar_stl_asthma
, which contains asthma rates for each census tract in the city. The interpolation process is very similar to the spatially extensive workflow, except with how the areal weight is calculated. Instead of using the source data’s area for reference, the target data is used in the denominator. Let:
Spatially intensive operations are used when the data to be interpolated are, for example, a percentage or density value. An example of these data can be found in ar_stl_asthma
, which contains asthma rates for each census tract in the city. The interpolation process is very similar to the spatially extensive workflow, except with how the areal weight is calculated. Instead of using the source data’s area for reference, the target data is used in the denominator. Let:
We can calculate the intensive interpolation by specifying a variable name for the intensive
argument in aw_interpolate()
and omitting the extensive
argument:
aw_interpolate(wards, tid = WARD, source = asthma, sid = GEOID,
weight = "sum", output = "tibble", intensive = "ASTHMA")
-#> # A tibble: 28 x 2
-#> WARD ASTHMA
-#> <int> <dbl>
-#> 1 1 13.4
-#> 2 2 13.2
-#> 3 3 14.1
-#> 4 4 13.6
-#> 5 5 13.8
-#> 6 6 11.7
-#> 7 7 9.72
-#> 8 8 9.82
-#> 9 9 11.8
-#> 10 10 9.44
+#> # A tibble: 28 x 4
+#> OBJECTID WARD AREA ASTHMA
+#> <dbl> <int> <dbl> <dbl>
+#> 1 1 1 46138761. 13.4
+#> 2 2 2 267817711. 13.2
+#> 3 3 3 66291644. 14.1
+#> 4 4 4 53210707. 13.6
+#> 5 5 5 60462396. 13.8
+#> 6 6 6 64337271. 11.7
+#> 7 7 7 101268146. 9.72
+#> 8 8 8 45966410. 9.82
+#> 9 9 9 73993891. 11.8
+#> 10 10 10 62915358. 9.44
#> # … with 18 more rows
This gives us an estimate of the asthma rates at the ward level.
Users should take care to consider the implications of interpolating multiple values, such as total population and the number of African American residents (both extensive), and then calculating a spatially intensive variable from them such as percent African American. Doing so in multiple steps, and thereby treating extensive and intensive values as independent, may result in estimates that differ from a single step process where the percent of African American residents is interpoltated directly.
+# re-load data
+race <- ar_stl_race
+
+# create combined data
+race %>%
+ select(GEOID, WHITE_E, BLACK_E) %>%
+ mutate(
+ TOTAL = WHITE_E+BLACK_E,
+ WHITE_PCT = WHITE_E/TOTAL,
+ BLACK_PCT = BLACK_E/TOTAL,
+ TOTAL_PCT = WHITE_PCT+BLACK_PCT
+ ) -> constrainedData
+
+# interpolate
+result <- aw_interpolate(ar_stl_wards, tid = WARD,
+ source = constrainedData, sid = GEOID,
+ weight = "sum", output = "tibble",
+ intensive = c("WHITE_PCT", "BLACK_PCT", "TOTAL_PCT"),
+ extensive = c("TOTAL", "WHITE_E", "BLACK_E"))
+
+# calculate new percentages
+result %>%
+ mutate(
+ WHITE_PCT_2 = WHITE_E/TOTAL,
+ BLACK_PCT_2 = BLACK_E/TOTAL,
+ TOTAL_PCT_2 = WHITE_PCT_2+BLACK_PCT_2
+ ) -> result
+
+# display
+result %>%
+ select(WHITE_PCT, WHITE_PCT_2, BLACK_PCT, BLACK_PCT_2, TOTAL_PCT, TOTAL_PCT_2)
+#> # A tibble: 28 x 6
+#> WHITE_PCT WHITE_PCT_2 BLACK_PCT BLACK_PCT_2 TOTAL_PCT TOTAL_PCT_2
+#> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
+#> 1 0.0195 0.0193 0.980 0.981 1.000 1
+#> 2 0.177 0.111 0.823 0.889 1 1.000
+#> 3 0.104 0.0818 0.896 0.918 1 1
+#> 4 0.0252 0.0191 0.975 0.981 1 1
+#> 5 0.174 0.178 0.826 0.822 1 1.000
+#> 6 0.409 0.440 0.591 0.560 1 1.
+#> 7 0.575 0.554 0.425 0.446 1.000 1
+#> 8 0.700 0.667 0.300 0.333 1 1
+#> 9 0.505 0.516 0.495 0.484 1 1
+#> 10 0.862 0.839 0.138 0.161 1 1.
+#> # … with 18 more rows
Note that there are a number of points of departure between the data interpolated as intensive values (WHITE_PCT
, BLACK_PCT
) and those that were interpolated as count data (i.e. extensive values) and then converted to intensive variables (WHITE_PCT_2
and BLACK_PCT_2
).
All of the above examples have created a tibble for output, but areal
also supports the creation of sf
objects as well:
aw_interpolate(wards, tid = WARD, source = asthma, sid = GEOID,
- weight = "sum", output = "sf", intensive = "ASTHMA")
-#> Simple feature collection with 28 features and 4 fields
-#> geometry type: POLYGON
-#> dimension: XY
-#> bbox: xmin: 733361.8 ymin: 4268336 xmax: 746157.7 ymax: 4295504
-#> epsg (SRID): 26915
-#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
-#> First 10 features:
-#> OBJECTID WARD AREA ASTHMA geometry
-#> 1 1 1 46138761 13.379577 POLYGON ((740184.2 4286431,...
-#> 2 2 2 267817711 13.164879 POLYGON ((742392.1 4289178,...
-#> 3 3 3 66291644 14.108107 POLYGON ((742956.1 4284113,...
-#> 4 4 4 53210707 13.613989 POLYGON ((739557.6 4284080,...
-#> 5 5 5 60462396 13.789757 POLYGON ((744883.8 4281632,...
-#> 6 6 6 64337271 11.747642 POLYGON ((742501.6 4279976,...
-#> 7 7 7 101268146 9.722839 POLYGON ((745618.6 4279867,...
-#> 8 8 8 45966410 9.823627 POLYGON ((739842.8 4277724,...
-#> 9 9 9 73993891 11.820622 POLYGON ((742619.4 4276734,...
-#> 10 10 10 62915358 9.439929 POLYGON ((737257.7 4277050,...
aw_interpolate(wards, tid = WARD, source = asthma, sid = GEOID,
+ weight = "sum", output = "sf", intensive = "ASTHMA")
+#> Simple feature collection with 28 features and 4 fields
+#> geometry type: POLYGON
+#> dimension: XY
+#> bbox: xmin: 733361.8 ymin: 4268336 xmax: 746157.7 ymax: 4295504
+#> epsg (SRID): 26915
+#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
+#> First 10 features:
+#> OBJECTID WARD AREA ASTHMA geometry
+#> 1 1 1 46138761 13.379577 POLYGON ((740184.2 4286431,...
+#> 2 2 2 267817711 13.164879 POLYGON ((742392.1 4289178,...
+#> 3 3 3 66291644 14.108107 POLYGON ((742956.1 4284113,...
+#> 4 4 4 53210707 13.613989 POLYGON ((739557.6 4284080,...
+#> 5 5 5 60462396 13.789757 POLYGON ((744883.8 4281632,...
+#> 6 6 6 64337271 11.747642 POLYGON ((742501.6 4279976,...
+#> 7 7 7 101268146 9.722839 POLYGON ((745618.6 4279867,...
+#> 8 8 8 45966410 9.823627 POLYGON ((739842.8 4277724,...
+#> 9 9 9 73993891 11.820622 POLYGON ((742619.4 4276734,...
+#> 10 10 10 62915358 9.439929 POLYGON ((737257.7 4277050,...
The sf
option will include all of the variables that were included in the original target data. The aw_interpolate()
function is pipe-able, allowing for existing tidyverse workflows to be integrated into the interpolation process. For example, if we wanted to remove the OBJECTID
and AREA
columns because they are not needed, this can be accomplished easily with areal
and dplyr
:
wards %>%
- select(-OBJECTID, -AREA) %>%
- aw_interpolate(tid = WARD, source = asthma, sid = GEOID,
- weight = "sum", output = "tibble", intensive = "ASTHMA")
-#> # A tibble: 28 x 2
-#> WARD ASTHMA
-#> <int> <dbl>
-#> 1 1 13.4
-#> 2 2 13.2
-#> 3 3 14.1
-#> 4 4 13.6
-#> 5 5 13.8
-#> 6 6 11.7
-#> 7 7 9.72
-#> 8 8 9.82
-#> 9 9 11.8
-#> 10 10 9.44
-#> # … with 18 more rows
All of the areal
functions that are exported support non-standard evaluation, meaning that inputs can be either unquoted as they are above or quoted:
wards %>%
select(-OBJECTID, -AREA) %>%
- aw_interpolate(tid = "WARD", source = asthma, sid = "GEOID",
+ aw_interpolate(tid = WARD, source = asthma, sid = GEOID,
weight = "sum", output = "tibble", intensive = "ASTHMA")
#> # A tibble: 28 x 2
#> WARD ASTHMA
@@ -576,6 +604,25 @@
#> 9 9 11.8
#> 10 10 9.44
#> # … with 18 more rows
All of the areal
functions that are exported support non-standard evaluation, meaning that inputs can be either unquoted as they are above or quoted:
wards %>%
+ select(-OBJECTID, -AREA) %>%
+ aw_interpolate(tid = "WARD", source = asthma, sid = "GEOID",
+ weight = "sum", output = "tibble", intensive = "ASTHMA")
+#> # A tibble: 28 x 2
+#> WARD ASTHMA
+#> <int> <dbl>
+#> 1 1 13.4
+#> 2 2 13.2
+#> 3 3 14.1
+#> 4 4 13.6
+#> 5 5 13.8
+#> 6 6 11.7
+#> 7 7 9.72
+#> 8 8 9.82
+#> 9 9 11.8
+#> 10 10 9.44
+#> # … with 18 more rows
This functionality is not available for the intensive
and extensive
arguments at this time.
areal
purposely exports the sub-functions that are called by aw_interpolate()
so that the interpolation process is not a “black box” but rather can be recreated manually. This is envisioned as a diagnostic toolkit, with the final interpolations estimated using the simpler aw_interpolate()
function once any issues have been identified and ameliorated
First, we’ll prepare the data but retaining only the columns we are interested in from the source data using the select()
function from dplyr
:
We want to be careful to retain both a column with a value to be interpolated (total population in this case, TOTAL_E
) and a column with a unique identification number (GEOID
in this case).
As we noted above, the interpolation process begins with calculating the intersection between the source and target data. We use the function aw_intersect()
to accomplish this:
wards %>%
- aw_intersect(source = race, areaVar = "area") -> intersect
-
-intersect
-#> Simple feature collection with 287 features and 4 fields
-#> geometry type: MULTIPOLYGON
-#> dimension: XY
-#> bbox: xmin: 733361.8 ymin: 4268411 xmax: 746155.7 ymax: 4295504
-#> epsg (SRID): 26915
-#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
-#> First 10 features:
-#> GEOID TOTAL_E WARD area geometry
-#> 1 29510106500 2863 1 4.271176e+05 MULTIPOLYGON (((738110.9 42...
-#> 2 29510106700 3168 1 5.946820e+05 MULTIPOLYGON (((738335.2 42...
-#> 3 29510107500 1991 1 6.146361e+05 MULTIPOLYGON (((740111.6 42...
-#> 4 29510107600 1958 1 1.071347e+06 MULTIPOLYGON (((739386.1 42...
-#> 5 29510110100 2820 1 5.563665e+01 MULTIPOLYGON (((739542.5 42...
-#> 6 29510109600 3050 1 3.586828e+03 MULTIPOLYGON (((740899.6 42...
-#> 7 29510126900 4462 1 8.903340e+05 MULTIPOLYGON (((739115.9 42...
-#> 8 29510127000 1826 1 1.015451e+03 MULTIPOLYGON (((740936.1 42...
-#> 9 29510106400 1986 1 6.866924e+05 MULTIPOLYGON (((738500.9 42...
-#> 10 29510108100 3309 1 1.198667e+02 MULTIPOLYGON (((740163.9 42...
wards %>%
+ aw_intersect(source = race, areaVar = "area") -> intersect
+
+intersect
+#> Simple feature collection with 287 features and 4 fields
+#> geometry type: GEOMETRY
+#> dimension: XY
+#> bbox: xmin: 733361.8 ymin: 4268411 xmax: 746155.7 ymax: 4295504
+#> epsg (SRID): 26915
+#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
+#> First 10 features:
+#> GEOID TOTAL_E WARD area geometry
+#> 1 29510106500 2863 1 4.271176e+05 POLYGON ((738110.9 4283989,...
+#> 2 29510106700 3168 1 5.946820e+05 POLYGON ((738335.2 4283055,...
+#> 3 29510107500 1991 1 6.146361e+05 POLYGON ((740111.6 4286547,...
+#> 4 29510107600 1958 1 1.071347e+06 POLYGON ((739386.1 4285082,...
+#> 5 29510110100 2820 1 5.563665e+01 POLYGON ((739542.5 4283728,...
+#> 6 29510109600 3050 1 3.586828e+03 POLYGON ((740899.6 4285240,...
+#> 7 29510126900 4462 1 8.903340e+05 POLYGON ((739115.9 4285807,...
+#> 8 29510127000 1826 1 1.015451e+03 MULTIPOLYGON (((740936.1 42...
+#> 9 29510106400 1986 1 6.866924e+05 POLYGON ((738500.9 4284750,...
+#> 10 29510108100 3309 1 1.198667e+02 POLYGON ((740163.9 4286467,...
Note that aw_intersect()
automatically calculates the area of the intersected feature.
Next, we want to apply the total area of our source features to our data using aw_total()
. This will implement the correct areal weighting approach based on the type
and weight
arguments. We’ll use the "sum"
approach to areal weights here:
intersect %>%
- aw_total(source = race, id = GEOID, areaVar = "area", totalVar = "totalArea",
- type = "extensive", weight = "sum") -> intersect
-
-intersect
-#> Simple feature collection with 287 features and 5 fields
-#> geometry type: MULTIPOLYGON
-#> dimension: XY
-#> bbox: xmin: 733361.8 ymin: 4268411 xmax: 746155.7 ymax: 4295504
-#> epsg (SRID): 26915
-#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
-#> First 10 features:
-#> GEOID TOTAL_E WARD area totalArea
-#> 1 29510106500 2863 1 4.271176e+05 1064776.2
-#> 2 29510106700 3168 1 5.946820e+05 1239241.5
-#> 3 29510107500 1991 1 6.146361e+05 912325.0
-#> 4 29510107600 1958 1 1.071347e+06 1146634.4
-#> 5 29510110100 2820 1 5.563665e+01 994606.5
-#> 6 29510109600 3050 1 3.586828e+03 3157835.4
-#> 7 29510126900 4462 1 8.903340e+05 4334701.5
-#> 8 29510127000 1826 1 1.015451e+03 14726832.5
-#> 9 29510106400 1986 1 6.866924e+05 1107191.4
-#> 10 29510108100 3309 1 1.198667e+02 3216478.2
-#> geometry
-#> 1 MULTIPOLYGON (((738110.9 42...
-#> 2 MULTIPOLYGON (((738335.2 42...
-#> 3 MULTIPOLYGON (((740111.6 42...
-#> 4 MULTIPOLYGON (((739386.1 42...
-#> 5 MULTIPOLYGON (((739542.5 42...
-#> 6 MULTIPOLYGON (((740899.6 42...
-#> 7 MULTIPOLYGON (((739115.9 42...
-#> 8 MULTIPOLYGON (((740936.1 42...
-#> 9 MULTIPOLYGON (((738500.9 42...
-#> 10 MULTIPOLYGON (((740163.9 42...
Changing type
to "intensive"
would be necessary for spatially intensive interpolations. Likewise, changing weight
to "total"
is necessary if areas that lack overlap should not be allocated into the target features.
With the total weight in hand, we are ready to calculate the areal weight itself using aw_weight()
.
intersect %>%
- aw_weight(areaVar = "area", totalVar = "totalArea",
- areaWeight = "areaWeight") -> intersect
+ aw_total(source = race, id = GEOID, areaVar = "area", totalVar = "totalArea",
+ type = "extensive", weight = "sum") -> intersect
intersect
-#> Simple feature collection with 287 features and 6 fields
-#> geometry type: MULTIPOLYGON
+#> Simple feature collection with 287 features and 5 fields
+#> geometry type: GEOMETRY
#> dimension: XY
#> bbox: xmin: 733361.8 ymin: 4268411 xmax: 746155.7 ymax: 4295504
#> epsg (SRID): 26915
#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
#> First 10 features:
-#> GEOID TOTAL_E WARD area totalArea areaWeight
-#> 1 29510106500 2863 1 4.271176e+05 1064776.2 4.011337e-01
-#> 2 29510106700 3168 1 5.946820e+05 1239241.5 4.798758e-01
-#> 3 29510107500 1991 1 6.146361e+05 912325.0 6.737030e-01
-#> 4 29510107600 1958 1 1.071347e+06 1146634.4 9.343406e-01
-#> 5 29510110100 2820 1 5.563665e+01 994606.5 5.593835e-05
-#> 6 29510109600 3050 1 3.586828e+03 3157835.4 1.135850e-03
-#> 7 29510126900 4462 1 8.903340e+05 4334701.5 2.053968e-01
-#> 8 29510127000 1826 1 1.015451e+03 14726832.5 6.895245e-05
-#> 9 29510106400 1986 1 6.866924e+05 1107191.4 6.202111e-01
-#> 10 29510108100 3309 1 1.198667e+02 3216478.2 3.726645e-05
+#> GEOID TOTAL_E WARD area totalArea
+#> 1 29510106500 2863 1 4.271176e+05 1064776.2
+#> 2 29510106700 3168 1 5.946820e+05 1239241.5
+#> 3 29510107500 1991 1 6.146361e+05 912325.0
+#> 4 29510107600 1958 1 1.071347e+06 1146634.4
+#> 5 29510110100 2820 1 5.563665e+01 994606.5
+#> 6 29510109600 3050 1 3.586828e+03 3157835.4
+#> 7 29510126900 4462 1 8.903340e+05 4334701.5
+#> 8 29510127000 1826 1 1.015451e+03 14726832.5
+#> 9 29510106400 1986 1 6.866924e+05 1107191.4
+#> 10 29510108100 3309 1 1.198667e+02 3216478.2
#> geometry
-#> 1 MULTIPOLYGON (((738110.9 42...
-#> 2 MULTIPOLYGON (((738335.2 42...
-#> 3 MULTIPOLYGON (((740111.6 42...
-#> 4 MULTIPOLYGON (((739386.1 42...
-#> 5 MULTIPOLYGON (((739542.5 42...
-#> 6 MULTIPOLYGON (((740899.6 42...
-#> 7 MULTIPOLYGON (((739115.9 42...
+#> 1 POLYGON ((738110.9 4283989,...
+#> 2 POLYGON ((738335.2 4283055,...
+#> 3 POLYGON ((740111.6 4286547,...
+#> 4 POLYGON ((739386.1 4285082,...
+#> 5 POLYGON ((739542.5 4283728,...
+#> 6 POLYGON ((740899.6 4285240,...
+#> 7 POLYGON ((739115.9 4285807,...
#> 8 MULTIPOLYGON (((740936.1 42...
-#> 9 MULTIPOLYGON (((738500.9 42...
-#> 10 MULTIPOLYGON (((740163.9 42...
Changing type
to "intensive"
would be necessary for spatially intensive interpolations. Likewise, changing weight
to "total"
is necessary if areas that lack overlap should not be allocated into the target features.
We can then multiply the value (TOTAL_E
) by the weight (areaWeight
) to get a population estimate for each intersected feature using aw_calculate()
:
With the total weight in hand, we are ready to calculate the areal weight itself using aw_weight()
.
intersect %>%
- aw_calculate(value = TOTAL_E, areaWeight = "areaWeight") -> intersect
-
-intersect
-#> Simple feature collection with 287 features and 6 fields
-#> geometry type: MULTIPOLYGON
-#> dimension: XY
-#> bbox: xmin: 733361.8 ymin: 4268411 xmax: 746155.7 ymax: 4295504
-#> epsg (SRID): 26915
-#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
-#> First 10 features:
-#> GEOID TOTAL_E WARD area totalArea areaWeight
-#> 1 29510106500 1148.4457826 1 4.271176e+05 1064776.2 4.011337e-01
-#> 2 29510106700 1520.2464178 1 5.946820e+05 1239241.5 4.798758e-01
-#> 3 29510107500 1341.3426720 1 6.146361e+05 912325.0 6.737030e-01
-#> 4 29510107600 1829.4388297 1 1.071347e+06 1146634.4 9.343406e-01
-#> 5 29510110100 0.1577461 1 5.563665e+01 994606.5 5.593835e-05
-#> 6 29510109600 3.4643427 1 3.586828e+03 3157835.4 1.135850e-03
-#> 7 29510126900 916.4807234 1 8.903340e+05 4334701.5 2.053968e-01
-#> 8 29510127000 0.1259072 1 1.015451e+03 14726832.5 6.895245e-05
-#> 9 29510106400 1231.7392274 1 6.866924e+05 1107191.4 6.202111e-01
-#> 10 29510108100 0.1233147 1 1.198667e+02 3216478.2 3.726645e-05
-#> geometry
-#> 1 MULTIPOLYGON (((738110.9 42...
-#> 2 MULTIPOLYGON (((738335.2 42...
-#> 3 MULTIPOLYGON (((740111.6 42...
-#> 4 MULTIPOLYGON (((739386.1 42...
-#> 5 MULTIPOLYGON (((739542.5 42...
-#> 6 MULTIPOLYGON (((740899.6 42...
-#> 7 MULTIPOLYGON (((739115.9 42...
-#> 8 MULTIPOLYGON (((740936.1 42...
-#> 9 MULTIPOLYGON (((738500.9 42...
-#> 10 MULTIPOLYGON (((740163.9 42...
There is an optional newVar
argument that can be used to store the estimates in a new column rather than in the existing value
column.
Finally, we aggregate the estimated values by target features using aw_aggregate()
:
We can then multiply the value (TOTAL_E
) by the weight (areaWeight
) to get a population estimate for each intersected feature using aw_calculate()
:
intersect %>%
- aw_aggregate(target = wards, tid = WARD, interVar = TOTAL_E) -> result
+ aw_calculate(value = TOTAL_E, areaWeight = "areaWeight") -> intersect
-result
-#> Simple feature collection with 28 features and 2 fields
-#> geometry type: POLYGON
+intersect
+#> Simple feature collection with 287 features and 6 fields
+#> geometry type: GEOMETRY
#> dimension: XY
-#> bbox: xmin: 733361.8 ymin: 4268336 xmax: 746157.7 ymax: 4295504
+#> bbox: xmin: 733361.8 ymin: 4268411 xmax: 746155.7 ymax: 4295504
#> epsg (SRID): 26915
#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
#> First 10 features:
-#> WARD TOTAL_E geometry
-#> 1 1 7991.565 POLYGON ((740184.2 4286431,...
-#> 2 2 12145.021 POLYGON ((742392.1 4289178,...
-#> 3 3 7344.287 POLYGON ((742956.1 4284113,...
-#> 4 4 8457.672 POLYGON ((739557.6 4284080,...
-#> 5 5 8783.377 POLYGON ((744883.8 4281632,...
-#> 6 6 14050.399 POLYGON ((742501.6 4279976,...
-#> 7 7 15840.086 POLYGON ((745618.6 4279867,...
-#> 8 8 12188.131 POLYGON ((739842.8 4277724,...
-#> 9 9 14217.149 POLYGON ((742619.4 4276734,...
-#> 10 10 11239.213 POLYGON ((737257.7 4277050,...
There is an optional newVar
argument that can be used to store the estimates in a new column rather than in the existing value
column.
Finally, we aggregate the estimated values by target features using aw_aggregate()
:
intersect %>%
+ aw_aggregate(target = wards, tid = WARD, interVar = TOTAL_E) -> result
+
+result
+#> Simple feature collection with 28 features and 2 fields
+#> geometry type: POLYGON
+#> dimension: XY
+#> bbox: xmin: 733361.8 ymin: 4268336 xmax: 746157.7 ymax: 4295504
+#> epsg (SRID): 26915
+#> proj4string: +proj=utm +zone=15 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
+#> First 10 features:
+#> WARD TOTAL_E geometry
+#> 1 1 7991.565 POLYGON ((740184.2 4286431,...
+#> 2 2 12145.021 POLYGON ((742392.1 4289178,...
+#> 3 3 7344.287 POLYGON ((742956.1 4284113,...
+#> 4 4 8457.672 POLYGON ((739557.6 4284080,...
+#> 5 5 8783.377 POLYGON ((744883.8 4281632,...
+#> 6 6 14050.399 POLYGON ((742501.6 4279976,...
+#> 7 7 15840.086 POLYGON ((745618.6 4279867,...
+#> 8 8 12188.131 POLYGON ((739842.8 4277724,...
+#> 9 9 14217.149 POLYGON ((742619.4 4276734,...
+#> 10 10 11239.213 POLYGON ((737257.7 4277050,...
R
vignettes/areal.Rmd
areal.Rmd
vignettes/data-preparation.Rmd
data-preparation.Rmd
NEWS.md
+ Patch fix to Issue 6 so that data not requiring the geometry collection fix do not get processed through that code, improving performance for those interpolations.
Updated draft of JOSS paper manuscript along with appendix code and response to reviewers added in paper/
.