distRcpp/README.Rmd at main · btskinner/distRcpp · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
---
title: distRcpp
output: md_document
---

```{r, include = FALSE}
options(width = 100)
```

![](https://img.shields.io/badge/dev-beta-red.svg)
[![GitHub release](https://img.shields.io/github/release/btskinner/distRcpp.svg)](https://github.com/btskinner/distRcpp)
[![R build status](https://github.com/btskinner/distRcpp/actions/workflows/check-standard.yml/badge.svg)](https://github.com/btskinner/distRcpp/actions)

This package uses [Rcpp](http://www.rcpp.org) to quickly compute population/distance-weighted measures. Geodesic distances can be computed using either [Haversine](https://en.wikipedia.org/wiki/Haversine_formula) or [Vincenty](https://en.wikipedia.org/wiki/Vincenty%27s_formulae) formulas. The package also has functions to return raw distance measures. If you are able to [install Rcpp on your machine](https://github.com/RcppCore/Rcpp), you should be able to install this package and use these functions.

Install the latest development version from Github with

```{r, eval = FALSE}
devtools::install_github("btskinner/distRcpp")
```

**NB** This package is still in early beta stages. It does not have much in the way of error handling. Data must be pre-processed so that no missing (`NA`) values are given to the functions.

## Available functions

### Weighted measures

#### `dist_weighted_mean()`

Interpolate values for a vector of locations (**X**) that are the inverse-distance-weighted average of measures taken at surrounding locations (**Y**). For each point, *x*, nearby values of the measure taken at **Y** are weighted more heavily than those from locations that are farther away.

#### `popdist_weighted_mean()`

Interpolate values for a vector of locations (**X**) that are the population/inverse-distance-weighted average of measures taken at surrounding locations (**Y**). For each point, *x*, nearby values of the measure taken at **Y** are weighted more heavily than those from locations that are farther away. Measures taken in more heavily populated *y* are given more weight than those with lower populations. This weighting scheme is a compromise between distance and population and is useful for interpolating measures that need to take both into account.

### Distances

#### `dist_1to1()`

Compute and return the geodesic distance between two spatial points. Returns distance in meters.

#### `dist_1tom()`

Compute and return the geodesic distance between one location and a vector of other locations. Returns vector of distances in meters.

#### `dist_mtom()`

Compute and return the geodesic distance between each coordinate pair in two vectors. Returns *n x k* matrix of distances in meters, where *n* = # of locations in first vector and *k* = # of locations in second vector.

#### `dist_df()`

Compute distance between corresponding coordinate pairs and return vector of distances in meters. For use when creating a new `data.frame` or [dplyr](https://CRAN.R-project.org/package=dplyr) `tbl_df()` column.

#### `dist_min()`

Compute minimum distance between each starting point, *x*, and
possible end points, **Y**. Returns vector of minimum distances in
meters that equals # of starting points (size of **X**).

#### `dist_max()`

Compute maximum distance between each starting point, *x*, and possible end points, **Y**. Returns vector of maximum distances in meters that equals # of starting points (size of **X**).

## Benchmark

Compare speed with base R function when measuring the distance between every United States population-weighted county centroid as measured in 2010 (N = 3,143 with complete measurements).

### Load data
```{r, include = FALSE}
options(tidyverse.quiet = TRUE)
```

```{r, results = 'hide', message = FALSE}
## libraries
libs <- c("tidyverse","microbenchmark","geosphere","distRcpp")
sapply(libs, require, character.only = TRUE)
```
```{r}
## read data
df <- get(data(county_centers))
df

## subset to 2010 population-weighted centroids (pclon10, pclat10)
p <- df %>% select(pclon10, pclat10) %>% drop_na() %>% data.frame()

```

### Check for equality
```{r}
dist_R <- geosphere::distm(p, p, fun = distHaversine)
dist_Rcpp <- distRcpp::dist_mtom(p[,1],p[,2],p[,1],p[,2])

dist_R[1:5,1:5]
dist_Rcpp[1:5,1:5]
all.equal(dist_R, dist_Rcpp)
```

### Benchmark

2018 MacBookPro, 2.9 GHz Intel Core i9, 32 GB 2400 MHz DDR4 SDRAM

```{r}
microbenchmark(
    dist_R = geosphere::distm(p, p, fun = distHaversine),
    dist_Rcpp = distRcpp::dist_mtom(p[,1],p[,2],p[,1],p[,2]),
    times = 100
)
```

### Big file
```{r}
## get census block group centers of population
bg <- readr::read_csv("https://www2.census.gov/geo/docs/reference/cenpop2010/blkgrp/CenPop2010_Mean_BG.txt") %>%
    setNames(tolower(names(.))) %>%
    filter(statefp < 56) %>%
    mutate(id = paste0(statefp, countyfp, tractce, blkgrpce),
           lon = longitude,
           lat = latitude) %>%
    select(id, lon, lat) %>%
    drop_na()

ct <- get(data(county_centers)) %>%
    rename(id = fips,
           lon = pclon10,
           lat = pclat10) %>%
    drop_na()
bg
ct

system.time(dist_Rcpp <- distRcpp::dist_min(x_df = ct, y_df = bg))

```