itsleeds / od Goto Github PK

As discussed with @joeytalbot

FYI,
cartography::getLinkLayer() does something similar to stplanr::od2line() and od::od_to_sf() but is clearly slower...

library(od)
bench::mark(check = FALSE, max_iterations = 100,
            stplanr = stplanr::od2line(od_data_df, od_data_zones),
            cartography = cartography::getLinkLayer(x = od_data_zones,df = od_data_df),
            od = od_to_sfc(od_data_df, od_data_zones),
            od_sf1 = od_to_sf(od_data_df, od_data_zones),
            od_sf2 = od_to_sf(od_data_df, od_data_zones, package = "sf", crs = 4326)
)
#> # A tibble: 5 x 6
#>   expression       min   median `itr/sec` mem_alloc `gc/sec`
#>   <bch:expr>  <bch:tm> <bch:tm>     <dbl> <bch:byt>    <dbl>
#> 1 stplanr        7.2ms   7.41ms      134.    75.7MB     2.84
#> 2 cartography   7.95ms   8.09ms      121.   490.4KB     2.05
#> 3 od            2.69ms   2.74ms      362.   515.9KB     0   
#> 4 od_sf1        3.16ms   3.21ms      308.    25.5KB     3.11
#> 5 od_sf2        3.35ms   3.42ms      289.   104.9KB     2.92

^{Created on 2020-02-12 by the reprex package (v0.3.0)}

Currently od_to_sf() does not copy across the CRS system and this has to be respecified

Created on my account as prototype but from ropensci/stplanr#367 it seems there is substantial buy-in and it seems moving it to https://github.com/ITSLeeds may be the right way to go.

Making it part of a team should make the package easier to maintain, more prominent and more future-proof is my thinking. Thoughts @mem48, @karlropkins (are you even on our GH team?), @HarukoNakao and others in ITS Leeds?

See here: https://win-builder.r-project.org/Q79aW1rb3yrs/00check.log

Only minor issue is:

* checking R code for possible problems ... [8s] NOTE
Found if() conditions comparing class() to string:
File 'od/R/oneway.R': if (class(x) == "factor") ...
File 'od/R/oneway.R': if (class(y) == "factor") ...
File 'od/R/oneway.R': if (class(x) == "factor") ...
File 'od/R/oneway.R': if (class(y) == "factor") ...
Use inherits() (or maybe is()) instead.

Including bidirectional flows and intrazonal flows.

See Robinlovelace/simodels#3

Currently od_to_sf() allows lines where one of the endpoints has NA coordinates

This is a sub-task of #20, with a view to splitting the task into small and manageable sub-tasks. You up for giving it a bash @mtennekes?

@pachamaltese is the maintainer of the really useful gravity package which provides a perfect partner for this package idea - this package could then provide a perfect interface and data-representation system to directly and simply fit a variety of interaction models onto OD data. gravity current accepts arbitrary input and requires users to specify which columns represent the flow and distance columns, along with possible additional regressors, so interface would be easily extended to here.

With #38 od_disaggregate() is more efficient, and deserves to be used in od_jitter().

The opposite of od_split().

Heads-up @agila5 can you paste a reproducible example illustrating the issue?

Hi Robin (cc @agila5 and @rogerbeecham),

Great to see your progress on the od package!

Finally, I have made some progress on 'my' od package, which I called odf. As you will see below, odf doen't have as much features as od. I have put the emphasis on the data structure, for which I used a graph structure. However, I have a different approach than I had in mind during my visit in Leeds: now I'm following a much more light-weighted approach.

Could you take a look and see if it's useful for your package? I've rendered a pdf from the vigette and also uploaded our draft paper. See also the latest version of the 'halfline' doughnut visualization, for which I've partly used odf.

I have the feeling that we can merge our approaches quite well. Let me know what you think.

Example from @joeytalbot: calculating desire lines for all DZs in Scotland.

Use case from @mem48: generate desire lines from LSOAs in England (32000 x 32000 = 1 billion rows).

E.g. from project with @joeytalbot 👍

od_intrazonal_disag1 = od::od_disaggregate(od = od_intrazonal_original, z = zones, subpoints = highways_osm_centroids_ed, population_column = "all", population_per_od = 100)
although coordinates are longitude/latitude, st_intersects assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
33 locations outside zones removed
although coordinates are longitude/latitude, st_intersects assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
Error in floor(x) : non-numeric argument to mathematical function
In addition: Warning message:
In Ops.factor(x, nrow(od_new)) : ‘/’ not meaningful for factors

Heads-up @natesheehan I'm working on this as fix for the bg traffic issue. Reproducible code:

# Aim: test the od_disaggregate() function

library(tidyverse)
library(tmap)

u = "https://github.com/cyipt/actdev/raw/main/data-small/poundbury/desire-lines-many.geojson"
desire_lines_many = sf::read_sf(u)

zones_msoa_national = pct::get_pct(layer = "z", geography = "msoa", national = TRUE)
centroids_lsoa_national = pct::get_pct(layer = "c", national = TRUE)
zones_lsoa_national = pct::get_pct(layer = "z", national = TRUE, geography = "lsoa")
zones_many = zones_msoa_national %>% filter(geo_code %in% desire_lines_many$geo_code2)
centroids_lsoa_many = centroids_lsoa_national[zones_many, ]
zones_lsoa_many = zones_lsoa_national %>% filter(geo_code %in% centroids_lsoa_many$geo_code)

qtm(zones_many, borders.lwd = 3) +
  qtm(desire_lines_many) +
  tm_shape(zones_lsoa_many) + tm_borders("red", lty = 3)

desire_lines_many = od::od_disaggregate(od = desire_lines_many %>% select(geo_code1:all_base), z = zones_many, subzones = zones_lsoa_many)

Result

See ropensci/stplanr#544

Happy to contribute this if you think it would be useful. The ideas is that you have points and want to make an OD datasets but cap the maximum distance considered. e.g. how we capped the PCT at 30 km.

Would need to add geodists as a dependancy.

Artist's impression:

With an argument to set maximum distances between points also, as discussed with @mem48

Issues / note to self

in Stplanr onewayid family can only sum. I've come accross a usecase for wanting another function.

I have some od data where the total flow between A and B is reported for both A > B and B > A, but due to counting errors the number are slightly different. In this case the max or min value would be needed not the sum, which would double the overall flow.

Also is na.rm supported?

library(od)
od = od_data_df[1:2, ]
zones = od::od_data_zones_min
subzones = od_data_zones_small
od_disag = od_disaggregate(od, zones, subzones)
#> Converting subzones to centroids
#> although coordinates are longitude/latitude, st_intersects assumes that they are planar
#> although coordinates are longitude/latitude, st_intersects assumes that they are planar
#> although coordinates are longitude/latitude, st_intersects assumes that they are planar
#> although coordinates are longitude/latitude, st_intersects assumes that they are planar
ncol(od_disag) -1 == ncol(od) # same number of columns (except disag data gained geometry)
#> [1] FALSE
sum(od_disag[[3]]) == sum(od[[3]])
#> [1] TRUE
sum(od_disag[[4]]) == sum(od[[4]])
#> [1] TRUE
plot(rowSums(sf::st_drop_geometry(od_disag)[4:10]), od_disag[[3]])

^{Created on 2021-04-30 by the reprex package (v2.0.0)}

when using od_aggregate, st_join sometimes joins 2 polygons to an o or d, due to the fact that aggzones are overlapping, and the o or d resides in both.
the function should throw an error when aggzones are intersecting.

I'll supply a reprex if requested

Thinking these are pretty fundamental functions for working with od data, for extracting only inter-zone and intra-zone pairs respectively. What do you think of this implementation @mem48 ?:

# todo: export at some point
od_interzone = function(x) {
 x[!x[[1]] == x[[2]], ]
}
od_intrazone = function(x) {
 x[x[[1]] == x[[2]], ]
}
library(od)
od_data = points_to_od(od_data_centroids)
nrow(od_data)
#> [1] 11449
nrow(od_interzone(od_data))
#> [1] 11342
nrow(od_intrazone(od_data))
#> [1] 107

^{Created on 2020-03-09 by the reprex package (v0.3.0)}

In the first instance points sampled from a route network but in the future simple random and sampled from polygons could be added (building footprints could be sampled from if the object using the proposed method).

Most importantly in the main od vignette https://itsleeds.github.io/od/articles/od.html and the README.

Example highlighting the issue:

set.seed(2021)
disag = od::od_disaggregate(od::od_data_df[1:2, ], z = od::od_data_zones_min)
#> Creating randomly sampled origin and destination points.
coords_origin = lwgeom::st_startpoint(disag)
#> Linking to GEOS 3.9.1, GDAL 3.3.2, PROJ 7.2.1
coords_destination = lwgeom::st_endpoint(disag)
summary(duplicated(coords_origin))
#>    Mode   FALSE    TRUE 
#> logical      26      17
summary(duplicated(coords_destination))
#>    Mode   FALSE    TRUE 
#> logical      29      14
plot(disag[0])

^{Created on 2021-11-28 by the reprex package (v2.0.1)}

This is documented here but should probably be in the intro vignette also.