Problemas de transformación de proyecciones poligonales dentro de la aplicación Shiny

Question

Recientemente he empezado a intentar crear aplicaciones de mapas con leaflet y shiny en R. Actualmente estoy perplejo por un problema con la proyección de polígonos dentro de shiny .

Aquí está el shapefile para reproducir mi código y aquí está el código despojado:

library(shiny)
library(leaflet) #MUST INSTALL DEVELOPMENT VERSION WITH:  devtools::install_github('rstudio/leaflet')
library(sp)
library(rgeos)
library(rgdal)

#pull in full rock country shapefile, set WGS84 CRS
countries <- readOGR("D:/NaturalEarth/HIF", layer = "ctry_hif", 
                     stringsAsFactors = F, encoding = "UTF-8")
countries <- spTransform(countries, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

#run shiny app
shinyApp(
  ui = fluidPage(leafletOutput('myMap', width = "80%", height = 500),
                 absolutePanel(width = "20%", top = 10, right = 5, 
                               selectInput(inputId = "location",
                                           label = "Country", 
                                           choices = c("", countries@data$sovereignt), 
                                           selected = "")
                 )
  ),     #END UI

  server <- function(input, output, session) {

    #observe click event; zooms to administrative units for the country that is clicked on
    observeEvent(input$location, {

      #subset countries shapefile by dropdown selection country
      selected <- countries[countries@data$sovereignt == input$location,]

      # selected <- spTransform(selected, CRS("+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
      # selected <- gBuffer(selected, width = -.2)
      # selected <- spTransform(selected, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

      #plot leaflet output
      output$myMap <- renderLeaflet({
        leaflet() %>%
          addTiles() %>%
          addPolygons(data = selected)
      })    #END TOP LEAFLET OUTPUT
    }) #END OBSERVE EVENT 
  } #END SERVER
) #END SHINYAPP

Ejecutar este código tal cual funciona bien: el mapa se actualiza con el país que el usuario seleccione en el menú desplegable.

El problema se produce cuando intento cambiar la proyección del selected polígono dentro del Shiny observeEvent. El objetivo final es poder colocar el polígono a nivel de país sobre un polígono a nivel regional con el rgeos función over y actualizar el mapa con regiones en lugar de países. He logrado esto dentro de Shiny para crear una aplicación perfectamente funcional. Pero cuando intento utilizar ese mismo código dentro de esta aplicación, no hay nada que hacer. He reducido el problema a las siguientes líneas de código (que están comentadas arriba):

selected <- spTransform(selected, CRS("+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
selected <- gBuffer(selected, width = -.2)
selected <- spTransform(selected, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

Cuando descomento esas líneas de código, la aplicación se rompe y me da el siguiente error:

Advertencia: Error en [[: subíndice fuera de los límites Rastreo de pila (más interno primero): 74: slot 73: .spTransform_Polygons 72: spTransformación 71: spTransformación 70: spTransformación 69: spTransformación 68: observeEventHandler [#19] 4: 3: do.call 2: print.shiny.appobj 1: ERROR: [on_request_read] conexión restablecida por el compañero

Cuando saco el código de Shiny y lo ejecuto con una cadena de países en lugar de input$location funciona perfectamente bien. Y de nuevo, he conseguido que el mismo código exacto se ejecute en una aplicación diferente con un evento de clic. Pero el cambio de proyección rompe esta.

¿Alguien tiene alguna idea de por qué este spTransform ¿se rompe mi aplicación Shiny?

Answer 1

Sólo añado un if a su código para filtrar la primera selección automática cuando se inicia la aplicación: "" .

Nota 1 : He seleccionado una mayor width en la opción gBuffer para ver cómo ha cambiado.

#observe click event; zooms to administrative units for the country that is clicked on

observeEvent(input$location, {

      #subset countries shapefile by dropdown selection country

      selected <- countries[countries@data$sovereignt == input$location,]

      if(input$location != "") {

        selected <- spTransform(selected, CRS("+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
        selected <- gBuffer(selected, width = -100000)
        selected <- spTransform(selected, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

      }

      #plot leaflet output
      output$myMap <- renderLeaflet({
        leaflet() %>%
          addTiles() %>%
          addPolygons(data = selected)

      })    #END TOP LEAFLET OUTPUT

    })

Si quieres que al iniciar la aplicación todos los países tengan el gBuffer efecto aplicado puede utilizar este código:

#observe click event; zooms to administrative units for the country that is clicked on
observeEvent(input$location, {

  #subset countries shapefile by dropdown selection country

  if(input$location == "" ) {

    selected <- countries

  } else {

    selected <- countries[countries@data$sovereignt == input$location,]

  }

  selected <- spTransform(selected, CRS("+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
  selected <- gBuffer(selected, width = -10000)
  selected <- spTransform(selected, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

  #plot leaflet output
  output$myMap <- renderLeaflet({
    leaflet() %>%
      addTiles() %>%
      addPolygons(data = selected)

  })    #END TOP LEAFLET OUTPUT

}) #END OBSERVE EVENT

Nota 2: Al principio no necesitas transformar tu objeto espacial porque tiene la proyección que deseabas desde el ctry_hif.prj archivo. Véase más abajo:

#pull in full rock country shapefile, set WGS84 CRS
countries <- readOGR(dsn = "D:/NaturalEarth/HIF", layer = "ctry_hif", stringsAsFactors = F, encoding = "UTF-8")

proj4string(countries)

# Output:
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"

# countries <- spTransform(countries, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")) # You don't need this