[R-sig-Geo] Match coordinates with NUTS3 regions

[Michele Molteni]

Hi everyone,

I have a .dta file with coordinates (WGS 84) and I would like to match them
to NUTS3 regions. I followed indications from a previous post but I get no
match even though the coordinates fall within NUTS3 bboxes. Could you
please be so kind to help me understand what I am doing wrong? This is what
I have so far:

> library(rgdal)
>
> ## Download Administrative Level data from EuroStat
>
> temp <- tempfile(fileext = ".zip")
> download.file("
https://ec.europa.eu/eurostat/cache/GISCO/distribution/v2/nuts/download/ref-nuts-2013-60m.shp.zip",
temp)
trying URL '
https://ec.europa.eu/eurostat/cache/GISCO/distribution/v2/nuts/download/ref-nuts-2013-60m.shp.zip
'
Content type 'application/zip' length 2330060 bytes (2.2 MB)
downloaded 2.2 MB
>
> unzip(temp)
>
> ##Read data after having unzipped it
> EU_NUTS <- readOGR(dsn = "C:/Users/User/Documents/TED_project/Replicating
paper/Shapefile", layer="NUTS_RG_60M_2013_4326_LEVL_3")
OGR data source with driver: ESRI Shapefile
Source: "C:\Users\User\Documents\TED_project\Replicating paper\Shapefile",
layer: "NUTS_RG_60M_2013_4326_LEVL_3"
with 1480 features
It has 5 fields
>
> ##Consider only EU29 countries (1361 NUTS3)
> EU_NUTS29 <- EU_NUTS[which(EU_NUTS using data$CNTR_CODE %in% c('BE', 'BG',
'CZ', 'DK', 'DE', 'EE', 'IE', 'EL', 'ES', 'FR', 'HR', 'IT', 'CY', 'LV',
'LT', 'LU', 'HU', 'MT', 'NL', 'AT', 'PL', 'PT', 'RO', 'SI', 'SK', 'FI',
'SE', 'UK', 'NO')), ]
>
> ##Import data contaning lat and lon and consider first 2500 observations
> awards_with_coordinates_first2500 <- read_dta("~/TED_project/Replicating
paper/awards_with_coordinates_first2500.dta")
> awards_with_coordinates_first2500 <- awards_with_coordinates_first2
500[(1:2500),]
>
> ## Here is the projection for map_nuts data
> proj4string(EU_NUTS29)
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
>
> ##Transform g_lat and g_lon in numeric
> awards_with_coordinates_first2500$g_lat <- as.numeric(as.character(awards
_with_coordinates_first2500$g_lat))
> awards_with_coordinates_first2500$g_lon <- as.numeric(as.character(awards
_with_coordinates_first2500$g_lon))
>
> ## Make a spatial points dataframe from point_data
> ## Assuming the projection to be longlat on the WGS84 datum
> ## First define the coordinates
> coordinates(awards_with_coordinates_first2500) <- ~g_lat + g_lon
>
> ## Then the projection
> proj4string(awards_with_coordinates_first2500) <- CRS("+proj=longlat
+datum=WGS84")
>
> ## Now transform the projection of awards_with_coordinates_first2500 to
that of EU_NUTS29
> awards_with_coordinates_first2500 <- spTransform(awards_with_coordinates_first2500,
proj4string(EU_NUTS29))
>
> ## Use the over() function to get the regions for the points
> over(awards_with_coordinates_first2500,EU_NUTS29,by.id=)
    NUTS_ID LEVL_CODE CNTR_CODE NUTS_NAME  FID
1      <NA>        NA      <NA>      <NA> <NA>
2      <NA>        NA      <NA>      <NA> <NA>

200    <NA>        NA      <NA>      <NA> <NA>
 [ reached 'max' / getOption("max.print") -- omitted 2300 rows ]
>
> bbox(awards_with_coordinates_first2500)
            min      max
g_lat 34.794891 51.46780
g_lon -3.018564 33.63785
> bbox(EU_NUTS29)
      min    max
x -61.806 55.826
y -21.350 71.127
> ##The points seem to be all in Africa
> plot(EU_NUTS29)
> plot(awards_with_coordinates_first2500, add = TRUE, col = "red", pch = 19)

Thanks in advance.

Cheers,
Michele.

	[[alternative HTML version deleted]]