Last active
March 20, 2022 23:46
-
-
Save mbacou/5880859 to your computer and use it in GitHub Desktop.
Working with attributes inside SpatialPolygonsDataFrame objects
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
## Example procedure for working with data (attributes) inside SpatialPolygonsDataFrame objects | |
# Load my favorite libraries for that sort of work | |
library(data.table) | |
library(rgdal) | |
# Load 2 shapefiles that, say, we want to merge. | |
# Note that you can read pretty much any spatial format using readOGR(). | |
tza.l2 <- readOGR("./analysis/TZA-AC-07/maps", "tza-ac-07_L2") | |
tza.l1 <- readOGR("./analysis/TZA-AC-07/maps", "tza-ac-07_L1") | |
# Let's explore the 1st object (for information) | |
str(tza.l2) | |
# Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots | |
# ..@ data :'data.frame': 140 obs. of 12 variables: | |
# .. ..$ ISO : Factor w/ 1 level "TZA": 1 1 1 1 1 1 1 1 1 1 ... | |
# .. ..$ NAME_0 : Factor w/ 1 level "Tanzania": 1 1 1 1 1 1 1 1 1 1 ... | |
# .. ..$ ID_1 : int [1:140] 3052 3052 3052 3052 3053 3053 3053 3054 3054 3054 ... | |
# .. ..$ NAME_1 : Factor w/ 26 levels "Arusha","Dar-Es-Salaam",..: 1 1 1 1 2 2 2 3 3 3 ... | |
# .. ..$ ID_2 : int [1:140] 33658 33659 33661 33662 33664 33665 33666 33667 33668 33669 ... | |
# .. ..$ NAME_2 : Factor w/ 130 levels "Arumeru","Arusha",..: 1 2 76 96 20 40 123 13 14 44 ... | |
# .. ..$ Shape_Leng: num [1:140] 2.801 0.504 7.466 6.13 1.063 ... | |
# .. ..$ Shape_Area: num [1:140] 0.2321 0.00882 1.30787 1.25441 0.02747 ... | |
# .. ..$ svyL1Code : int [1:140] 2 2 2 2 7 7 7 1 1 1 ... | |
# .. ..$ svyL1Name : Factor w/ 26 levels "Arusha","Dar es salaam",..: 1 1 1 1 2 2 2 3 3 3 ... | |
# .. ..$ svyL2Code : int [1:140] 7 3 1 5 2 1 3 6 5 1 ... | |
# .. ..$ svyL2Name : Factor w/ 136 levels "Arusha Rural",..: 1 2 81 101 22 44 130 5 16 48 ... | |
# ..@ polygons :List of 140 | |
# .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots | |
# .. .. .. ..@ Polygons :List of 49 | |
# .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots | |
# .. .. .. .. .. .. ..@ labpt : num [1:2] 36.91 -2.99 | |
# .. .. .. .. .. .. ..@ area : num 1.63e-16 | |
# .. .. .. .. .. .. ..@ hole : logi FALSE | |
# .. .. .. .. .. .. ..@ ringDir: int 1 | |
# .. .. .. .. .. .. ..@ coords : num [1:5, 1:2] 36.9 36.9 36.9 36.9 36.9 ... | |
# .. .. .. .. ..[list output truncated] | |
# .. .. .. ..@ plotOrder: int [1:49] 31 48 23 45 40 38 15 37 16 20 ... | |
# .. .. .. ..@ labpt : num [1:2] 36.73 -3.43 | |
# .. .. .. ..@ ID : chr "0" | |
# .. .. .. ..@ area : num 0.12 | |
# .. .. [list output truncated] | |
# ..@ plotOrder : int [1:140] 92 47 110 60 114 73 97 13 116 100 ... | |
# ..@ bbox : num [1:2, 1:2] 29.327 -11.746 40.445 -0.982 | |
# .. ..- attr(*, "dimnames")=List of 2 | |
# .. .. ..$ : chr [1:2] "x" "y" | |
# .. .. ..$ : chr [1:2] "min" "max" | |
# ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots | |
# .. .. ..@ projargs: chr NA | |
# That object is made of 5 "slots" | |
slotNames(tza.l2) | |
# [1] "data" "polygons" "plotOrder" "bbox" "proj4string" | |
# We can access each slot (using "@") | |
class(tza.l2@data) | |
# [1] "data.frame" | |
class(tza.l2@polygons) | |
# [1] "list" | |
# Each polygon is itself a list with 5 slots (incl. all the vertice coordinates). | |
slotNames(tza.l2@polygons[[1]]) | |
# [1] "Polygons" "plotOrder" "labpt" "ID" "area" | |
# The ID slot is what ties the geometries and the attributes. | |
all(lapply(tza.l2@polygons, slot, "ID") == row.names(tza.l2)) | |
# [1] TRUE | |
# Note that like row.names, polygon IDs are character not integer. | |
summary(as.integer(lapply(tza.l2@polygons, slot, "ID"))) | |
# Min. 1st Qu. Median Mean 3rd Qu. Max. | |
# 0.0 34.8 69.5 69.5 104.0 139.0 | |
summary(as.integer(row.names(tza.l2)) | |
# Min. 1st Qu. Median Mean 3rd Qu. Max. | |
# 0.0 34.8 69.5 69.5 104.0 139.0 | |
# Save the data row.names into an explicit variable | |
tza.l2$rn <- row.names(tza.l2) | |
# Create temporary data tables to work on the attributes | |
tmp.l2 <- data.table(tza.l2@data) | |
tmp.l1 <- data.table(tza.l1@data) | |
# Say we want to merge data from table tmp.l1 into tmp.l2 | |
setkey(tmp.l2, svyL1Name, svyL2Code) | |
setkey(tmp.l1, region, district) | |
tmp.l2 <- tmp.l1[tmp.l2] | |
# Then let's re-attach the table to the original SpatialPolygonsDataFrame | |
# (preserving the original order of the row.names) | |
setkey(tmp.l2, rn) | |
tza.l2@data <- tmp.l2[row.names(tza)] | |
# Export back to shapefile (or to any spatial format) | |
writeOGR(tza.l2, "./analysis/TZA-AC-07/maps", "tza-ac-07_L2", driver="ESRI Shapefile") |
I love you so mutch!!!!!!! Thank'ss
This is very useful.
I want to assign a value from a variable " VAR1"
, contained in another dataset that also contains the variable " ID_2"
, to each value of the variable "ID_2"
inside SpatialPolygonsDataFrame objects.
To accomplish this I replicated the code above except that I used inner_join
with by = "ID_2"
to add the new data (attribute) to the data of the SpatialPolygonsDataFrame object.
But when I re-attach the table to the original SpatialPolygonsDataFrame I get Error in setkeyv (x, cols, verbose = verbose, physical = physical): x is not a data.table.
I would appreciate any ideas.
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Thanks! Quick question, how would you go about subsetting any one of the polygons inside the object? For example, eliminating ID number 1 (TZA)?