R中QGIS工作时如何使用sf包中的read_sf()函数读取shape文件

Question

我有一个名为

gadm41_EGY_1.shp

的形状文件，我已成功使用它在 QGIS 中进行分析，我想使用

read_sf()

包中的函数

sf

将其读入 R。

我想在R中进行伪缺席分析，我需要读取形状文件并在R Studio的绘图环境中直观地绘制地图。

为了确保shape文件的直接路径正确，我手动设置了工作目录以确保路径正确。

有人知道我做错了什么吗？

1. 第一次尝试：

   #Install the pacakge sf
   install.packages("sf", configure.args = "--with-proj-lib=/usr/local/lib/")
   
   library(sf)
   
   #Set the CRS
   crswgs84=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
   
   #Plot the shape file 
   Egypt_Sf<-read_sf("/Users/Documents/gadm41_EGY_1_Shape/gadm41_EGY_1.shp", 
                       proj4string=crswgs84, verbose=TRUE)
   

   错误信息：

   Error: Cannot open "/Users/Documents/gadm41_EGY_1_Shape/gadm41_EGY_1.shp""; The source could be corrupt or not supported. See `st_drivers()` for a list of supported formats.
   In addition: Warning message:
   In CPL_read_ogr(dsn, layer, query, as.character(options), quiet,  :
     GDAL Error 4: Unable to open /Users/Documents/gadm41_EGY_1_Shape/gadm41_EGY_1.shp" or /Users/Documents/gadm41_EGY_1_Shape/gadm41_EGY_1.shp".SHX. Set SHAPE_RESTORE_SHX config option to YES to restore or create it.
   

2. 第二次尝试：

   #Define a path to the shapefile 
   Ds_sf <- system.file("/Users/Documents/gadm41_EGY_1_Shape/gadm41_EGY_1.shp")
   Ds_sf
   
   #Set the CRS
   crswgs84=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
   
   #Read shape file in R
   Egypt_Sf<-read_sf(Ds_sf, proj4string=crswgs84, quiet = TRUE, verbose=TRUE)
   

   错误信息：

   Error: `dsn` must point to a source, not an empty string.
   

3. 第三次尝试：

   Egypt_Sf <- read_sf(dsn = "/Users/Documents/gadm41_EGY_1_Shape/gadm41_EGY_1.shp", 
                            quiet=TRUE, verbose=TRUE, stringAsFactors=TRUE)
   

错误信息：

Error in st_sf(x, ..., agr = agr, sf_column_name = sf_column_name) : 
  no simple features geometry column present

Answer 1

您似乎正在处理来自 https://gadm.org/data.html 的 Shapefile，因为我们真的不知道您下载和处理文件的情况如何，从干净的副本开始可能更容易。

将 Shapefile 存档保存为

.shp.zip

时，

sf

/ GDAL 可以猜测它是压缩的 ESRI shapefile，如果可能包含多个数据集（如 GADM 数据集一样），并且这些数据集可以作为图层处理：

library(sf)
#> Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUE

gadm_egy_url <- "https://geodata.ucdavis.edu/gadm/gadm4.1/shp/gadm41_EGY_shp.zip"

# download, save as .shp.zip
download.file(gadm_egy_url, "gadm41_EGY.shp.zip")

# list layers (individual shp files in the archive)
st_layers("gadm41_EGY.shp.zip")
#> Driver: ESRI Shapefile 
#> Available layers:
#>     layer_name geometry_type features fields crs_name
#> 1 gadm41_EGY_0       Polygon        1      2   WGS 84
#> 2 gadm41_EGY_1       Polygon       27     11   WGS 84
#> 3 gadm41_EGY_2       Polygon      343     13   WGS 84

# load one 
read_sf("gadm41_EGY.shp.zip", layer = "gadm41_EGY_1")
#> Simple feature collection with 27 features and 11 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: 24.6981 ymin: 21.72539 xmax: 36.24875 ymax: 31.66792
#> Geodetic CRS:  WGS 84
#> # A tibble: 27 × 12
#>    GID_1  GID_0 COUNTRY NAME_1 VARNAME_1 NL_NAME_1 TYPE_1 ENGTYPE_1 CC_1  HASC_1
#>    <chr>  <chr> <chr>   <chr>  <chr>     <chr>     <chr>  <chr>     <chr> <chr> 
#>  1 EGY.1… EGY   Egypt   Ad Da… Al Daqah… الدقهلية  Muhaf… Governor… NA    EG.DQ 
#>  2 EGY.2… EGY   Egypt   Al Ba… Mar Rojo… محافظة ا… Muhaf… Governor… NA    EG.BA 
#>  3 EGY.3… EGY   Egypt   Al Bu… Beheira|… البحيرة   Muhaf… Governor… NA    EG.BH 
#>  ...

通过

sf

，我们还可以使用 GDAL 虚拟文件系统，例如定义我们正在处理 zip 文件的 URL，并让 GDAL 处理下载和解压：

read_sf(paste0("/vsizip/vsicurl/", gadm_egy_url, "/gadm41_EGY_1.shp"))
#> Simple feature collection with 27 features and 11 fields
#> Geometry type: MULTIPOLYGON
#> ... 

或者直接使用

geodata

pacakge 下载 GADM 数据集：

geodata::gadm("Egypt", path = "./") |> 
  st_as_sf()
#> Simple feature collection with 27 features and 11 fields
#> Geometry type: GEOMETRY
#> ...

^{创建于 2024 年 10 月 5 日，使用 reprex v2.1.1}