sf::st_read_db() speed compared to postGIStools::get_postgis_query()

[andschar]

I have a large (extent not size) polygon layer of the German federal state borders on my PostgreSQL data base. There is also an index on the geom column.

CREATE INDEX germany_bld_geom_idx
  ON spatial_derived.germany_bld
  USING gist
  (geom);

Using the sf:: package to load the layer into R takes extremely long ~40s:

# bld = sf::st_read_db(con, table = c('spatial_derived', 'germany_bld'))
bld = sf::st_read_db(con, query = "SELECT * FROM spatial_derived.germany_bld", geom_column = 'geom')

Using the postGIStools:: package, it is done in ~0.5s:

bld = postGIStools::get_postgis_query(con, "SELECT * FROM spatial_derived.germany_bld", geom_name="geom")
bld = sf::st_as_sf(bld)

Up to now I've only seen advantages in using sf:: in R, however in this case other packages are faster. Is st_read_db() still in development? What causes this big difference?

PS: The question is also on: https://gis.stackexchange.com/questions/257576/sfst-read-db-considerably-slower-than-postgistoolsget-postgis-query?noredirect=1#comment406825_257576

[edzer]

@etiennebr ? Since you're asking for the whole table, I don't think the index comes into play. Both routines seem to use dbGetQuery to query the database. Can you provide a reproducible example?

[andschar]

reproducible example with a DB running will be difficult but here's a link to the shape file I imported into the DB: http://www.geodatenzentrum.de/auftrag1/archiv/vektor/vg2500/2011/ (first one!; unziped: ./vg2500_2011-01-01.geo84.shape/vg2500/vg2500_bld.shp) and a R script how I uploaded and imported the layer:

require(sf)
require(RPostgreSQL)
require(postGIStools)
source('/home/andreas/Downloads/DB_access.R') # your access information!

# read shape file
bld = st_read('/home/andreas/Downloads/vg2500_2011-01-01.geo84.shape/vg2500/vg2500_bld.shp', stringsAsFactors = FALSE)

# change characters
bld[7,3] = 'Baden-Wuerttemberg'
bld[15,3] = 'Thueringen'
names(bld) = tolower(names(bld))

# transform
st_crs(bld)
bld2 = st_transform(bld, 31467)

# upload to DB
drv = dbDriver("PostgreSQL")
con = dbConnect(drv, dbname = DBname_1, user = DBuser, host = DBhost, port = DBport, password = DBpassword)

dbSendQuery(con, "DROP TABLE IF EXISTS spatial_derived.germany_bld;")
st_write_db(con, bld2,
            table = c("spatial_derived", "germany_bld"),
            geom_name = "geom",
            row.names = FALSE)
dbSendQuery(con, "CREATE INDEX germany_bld_geom_idx
                  ON spatial_derived.germany_bld
                  USING gist
                  (geom);")
dbSendQuery(con, "COMMENT ON TABLE spatial_derived.germany_bld IS 'Germany Bundeslaender';")

dbDisconnect(con)
dbUnloadDriver(drv)

# load data from DB
drv = dbDriver("PostgreSQL")
con = dbConnect(drv, dbname = DBname_1, user = DBuser, host = DBhost, port = DBport, password = DBpassword)

time = Sys.time()
bld_db1 = st_read_db(con, query = "SELECT * FROM spatial_derived.germany_bld",
                 geom_column = 'geom')
Sys.time() - time # ~40s

time = Sys.time()
bld_db2 = get_postgis_query(con, "SELECT * FROM spatial_derived.germany_bld",
                        geom_name = 'geom')
bld_db2 = st_as_sf(bld_db2)
Sys.time() - time # 0.2 s

dbDisconnect(con)
dbUnloadDriver(drv)

[edzer]

Thanks! I'm seeing quite the opposite:

> time = Sys.time()
> bld_db1 = st_read_db(con, query = "SELECT * FROM germany_bld",
+                  geom_column = 'geom')
> Sys.time() - time # ~40s
Time difference of 0.06120801 secs
> 
> time = Sys.time()
> bld_db2 = get_postgis_query(con, "SELECT * FROM germany_bld",
+                         geom_name = 'geom')
> bld_db2 = st_as_sf(bld_db2)
> Sys.time() - time # 0.2 s
Time difference of 0.2348094 secs

with the following modified script:

require(sf)
require(RPostgreSQL)
require(postGIStools)
#source('/home/andreas/Downloads/DB_access.R') # your access information!

# read shape file
bld = st_read('vg2500_2011-01-01.geo84.shape/vg2500/vg2500_bld.shp', stringsAsFactors = FALSE)

# change characters
bld[7,3] = 'Baden-Wuerttemberg'
bld[15,3] = 'Thueringen'
names(bld) = tolower(names(bld))

# transform
st_crs(bld)
bld2 = st_transform(bld, 31467)

# upload to DB
drv = dbDriver("PostgreSQL")
con = dbConnect(drv, dbname = "postgis")
#con = dbConnect(drv, dbname = DBname_1, user = DBuser, host = DBhost, port = DBport, password = DBpassword)

dbSendQuery(con, "DROP TABLE IF EXISTS germany_bld;")
st_write_db(con, bld2,
            table = c("germany_bld"),
            geom_name = "geom",
            row.names = FALSE)
dbSendQuery(con, "CREATE INDEX germany_bld_geom_idx
                  ON germany_bld
                  USING gist
                  (geom);")
dbSendQuery(con, "COMMENT ON TABLE germany_bld IS 'Germany Bundeslaender';")

dbDisconnect(con)
dbUnloadDriver(drv)

# load data from DB
drv = dbDriver("PostgreSQL")
con = dbConnect(drv, dbname = "postgis")

time = Sys.time()
bld_db1 = st_read_db(con, query = "SELECT * FROM germany_bld",
                 geom_column = 'geom')
Sys.time() - time # ~40s

time = Sys.time()
bld_db2 = get_postgis_query(con, "SELECT * FROM germany_bld",
                        geom_name = 'geom')
bld_db2 = st_as_sf(bld_db2)
Sys.time() - time # 0.2 s

dbDisconnect(con)
dbUnloadDriver(drv)

If you can explain me how to set up a schema in my database, I can try that as well.

[andschar]

you mean a schema in a PostgresSQL data base?
from R: dbSendQuery(con, "CREATE SCHEMA test")

[edzer]

Thanks, that didn't change: sf::st_read_db is still four times as fast as get_postgis_query.

[edzer]

sessionInfo()?

Mine:

> sessionInfo()
R version 3.4.2 (2017-09-28)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 16.04.3 LTS

Matrix products: default
BLAS: /usr/lib/openblas-base/libblas.so.3
LAPACK: /usr/lib/libopenblasp-r0.2.18.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_GB.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_GB.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] postGIStools_0.2.1 sp_1.2-5           RPostgreSQL_0.6-2  DBI_0.7           
[5] sf_0.5-5          

loaded via a namespace (and not attached):
 [1] compiler_3.4.2  magrittr_1.5    rgdal_1.2-12    tools_3.4.2    
 [5] units_0.4-6     Rcpp_0.12.13    udunits2_0.13   grid_3.4.2     
 [9] rgeos_0.3-25    lattice_0.20-35

[andschar]

Ok, I had some outdated packages. Updated them with CRAN and sessionInfo() looks like this now:

> sessionInfo()
R version 3.3.3 (2017-03-06)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Linux Mint 18

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=de_AT.UTF-8    LC_MESSAGES=en_US.UTF-8    LC_PAPER=de_AT.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=de_AT.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] postGIStools_0.2.1 sp_1.2-4           sf_0.5-4           RPostgreSQL_0.6-2  DBI_0.7           

loaded via a namespace (and not attached):
[1] magrittr_1.5    rgdal_1.2-12    tools_3.3.3     units_0.4-6     Rcpp_0.12.12    udunits2_0.13   grid_3.3.3     
[8] rgeos_0.3-25    lattice_0.20-34

I'm still at 40s. I guess the different R versions are not important here and all the geo and data base packages (postGIStools, sp, sf, RPostgreSQL, DBI, rgdal, rgeos) are at least up to date with the latest CRAN version. Did I overlook something?

[edzer]

Not only that!

[etiennebr]

I can't reproduce this. In my case, for instance, st_read_db() is consistantly 18 times faster than postGIStools::get_postgis_query() on a 8Gb linestring table on a remote server.

[dpprdan]

I also can't reproduce this. st_read_db is consistently 3-4 times faster than get_postgis_query, with similar timings as @edzer
PostgreSQL/PostGIS version or localhost/remote DB server does not matter, both have roughly the same timings.
R 3.4.1, sf 0.5.5,
localhost: PostgreSQL 9.6, PostGIS 2.3
remote: PostgreSQL 9.2, PostGIS 2.1.8