readOGR (rgdal) fails to fetch polygon names from

问题:

I am trying to import an KML map of CCG boundaries in England (Available here, 200Kb) into R using readOGR function from package rgdal. My end-goal is to create a heat-map by colouring CCGs according to some associated value. I have a list with those values next to CCG names in one data frame. I need to match CCG names in that data frame with CCG names in the imported map object, and assign colours based on the value. However, I cannot see any CCG names imported in the map object, although they are present in the KML file. This is what I am doing:

library(sp)
library(rgdal)
library(maps)
library(maptools)

Assuming the KML file is in the working directory. Listing layers:

ogrListLayers("Clinical_Commissioning_Groups_April_2016_Ultra_Generalised_Clipped_Boundaries_in_England.KML")

Reading OGRGeoJSON layer:

ccg_boundaries <- ReadOGR("Clinical_Commissioning_Groups_April_2016_Ultra_Generalised_Clipped_Boundaries_in_England.KML","OGRGeoJSON")

R Studio shows there are two sections (right word?) in the object.

polygons, which contains data for each polygon, e.g. for the first one:

> ccg_boundaries@polygons[1]
[[1]]
An object of class "Polygons"
Slot "Polygons":
[[1]]
An object of class "Polygon"
Slot "labpt":
[1] -2.104671 54.040320
Slot "area":
[1] 0.168067
...

And data, with two variables (Name and Description) which I would expect to contain CCG names, but it is empty:

> ccg_boundaries@data
    Name Description
0                   
1                   
2                   
3                   
4                   
5         

However, the CCG names are there in the KML file, which can be seen if opened with a Word editor, e.g. the first one in the alphabetic order is "NHS Airedale, Wharfedale and Craven".

<PolyStyle><fill>0</fill></PolyStyle></Style>
    <ExtendedData><SchemaData schemaUrl="#OGRGeoJSON">
        <SimpleData name="objectid">1</SimpleData>
        <SimpleData name="ccg16cd">E38000001</SimpleData>
        <SimpleData name="ccg16nm">NHS Airedale, Wharfedale and Craven CCG</SimpleData>

Is there maybe an option to readOGR or some other option to extract them and include in the object?

回答1:

OK, if anyone encounters the same problem, here is the solution I found.

The website provides the maps in two formats: KML and SHP. I chose KML, because this was used in a worked example that I was following. But there appears to be a problem with this particular KML file or how it was generated. I tried the procedure with a Shapefile (SHP) instead, and it worked like a charm.

Shapefiles can be read into R by the same function, but don't need specifying the layer:

ccg_boundaries <- ReadOGR("Clinical_Commissioning_Groups_April_2016_Ultra_Generalised_Clipped_Boundaries_in_England.SHP")

CCG names are now there in the ccg16nm variable:

> head(ccg_boundaries@data)
  objectid   ccg16cd                                 ccg16nm st_areasha st_lengths
0        1 E38000001 NHS Airedale, Wharfedale and Craven CCG 1224636590  193149.74
1        2 E38000002                         NHS Ashford CCG  582174805  122841.19
2        3 E38000003                  NHS Aylesbury Vale CCG  984352696  229544.11
3        4 E38000004            NHS Barking and Dagenham CCG   36315011   31196.87
4        5 E38000005                          NHS Barnet CCG   86654018   41833.69
5        6 E38000006                        NHS Barnsley CCG  327520495  106476.52

回答2:

Your issue is that windows does not have the necessary library to extract the ExtendedData from a KML.
I provided a working solution here: https://stackoverflow.com/a/51657844/2763996

The solution to your problem is the following function that will work on your example KML:

library(tidyverse)
library(xml2)
library(rgdal)

readKML <- function(file,keep_name_description=FALSE,layer,...) {
  # Set keep_name_description = TRUE to keep "Name" and "Description" columns
  #   in the resulting SpatialPolygonsDataFrame. Only works when there is
  #   ExtendedData in the kml file.

  sp_obj<-readOGR(file,layer,...)
  xml1<-read_xml(file)
  if (!missing(layer)) {
    different_layers <- xml_find_all(xml1, ".//d1:Folder") 
    layer_names <- different_layers %>% 
      xml_find_first(".//d1:name") %>% 
      xml_contents() %>% 
      xml_text()

    selected_layer <- layer_names==layer
    if (!any(selected_layer)) stop("Layer does not exist.")
    xml2 <- different_layers[selected_layer]
  } else {
    xml2 <- xml1
  }

  # extract name and type of variables

  variable_names1 <- 
    xml_find_first(xml2, ".//d1:ExtendedData") %>% 
    xml_children() 

  while(variable_names1 %>% 
        xml_attr("name") %>% 
        is.na() %>% 
        any()&variable_names1 %>%
        xml_children() %>% 
        length>0) variable_names1 <- variable_names1 %>%
    xml_children()

  variable_names <- variable_names1 %>%
    xml_attr("name") %>% 
    unique()

  # return sp_obj if no ExtendedData is present
  if (is.null(variable_names)) return(sp_obj)

  data1 <- xml_find_all(xml2, ".//d1:ExtendedData") %>% 
    xml_children()

  while(data1 %>%
        xml_children() %>% 
        length>0) data1 <- data1 %>%
    xml_children()

  data <- data1 %>% 
    xml_text() %>% 
    matrix(.,ncol=length(variable_names),byrow = TRUE) %>% 
    as.data.frame()

  colnames(data) <- variable_names

  if (keep_name_description) {
    sp_obj@data <- data
  } else {
    try(sp_obj@data <- cbind(sp_obj@data,data),silent=TRUE)
  }
  sp_obj
}