• Input: register data, demographic and mobility surveys
• Survey calibration: Compute weights
• Sampling without replacement: Transform to integer weights

• Statistical matching: Combine datasets

  

1. Weighted random sampling without replacement
2. Similarity-based statistical matching
3. Survey calibration

sample.int

function (n, size = n, replace = FALSE, prob = NULL) 
{
    if (!replace && is.null(prob) && n > 1e+07 && size <= n/2) 
        .Internal(sample2(n, size))
    else .Internal(sample(n, size, replace, prob))
}

• n: Number of items to choose from
• size: Number of items to choose
  • optimize for size ≫ 1
• replace: With or without replacement
• prob: Uniform or non-uniform probabilities

Common framework:

• Subdivide an interval according to probabilities
• Repeat
  • Uniformly sample point from interval
  • Choose sub-interval covered by point
  • If sampling without replacement, remove sub-interval

Common framework:

• Subdivide an interval according to probabilities
• Repeat
  • Uniformly sample point from interval
  • Choose sub-interval covered by point
  • If sampling without replacement, remove sub-interval

• R uses trivial algorithm with update in \(O(n)\)
• Heap-like data structure (Partial Sum Trees) with update in \(O(\log n)\) (Wong and Easton, 1980)
• Alternative approaches

• Rejection sampling
  • Sample with replacement, throw away duplicates
• One-pass sampling (Efraimidis and Spirakis, 2006)

wrswoR::sample_int_rank

function (n, size, prob) 
{
    .check_args(n, size, prob)
    head(order(rexp(n)/prob), size)
}

n %>% rexp %>% divide_by(prob) %>% order %>% head(size)

lm(log10(time) ~ a. + a.:log10(n) + a.:log10(size) - 1, weights = n)

                    Estimate Std. Error
a.crank               1.8789     0.0186
a.crank:log10(n)      0.7809     0.0043
a.crank:log10(size)   0.3106     0.0019
a.rej                 3.1159     0.0186
a.rej:log10(n)        0.4327     0.0043
a.rej:log10(size)     0.4296     0.0019
a.R                   0.2789     0.0186
a.R:log10(n)          1.0722     0.0043
a.R:log10(size)       0.8159     0.0019

10 ^ 0.02

[1] 1.047129

devtools::install_github("krlmlr/wrswoR")

• Input: two joint distributions \((X, Y)\) and \((X, Z)\)
• Output: a joint distribution \((X, Y, Z)\)

Hot deck: Combine observations from \((X, Y)\) with "matching" observations from \((X, Z)\) to create a realization of \((X, Y, Z)\).

\[(x, y) ⋈ (x', z) = (x, y, z)\]

• Precondition: \(Y\) and \(Z\) independent given \(X\)

• Variants:
  • Exact: \(x' = x\)
  • Approximate: \(x'\) similar to \(x\)

Distance metric for multivariate distributions (Gower, 1971)

Weighted sum of distances for each variable, in \([0, 1]\)

• Interval-scaled: Relative to interval width
• Ordinal variables: Relative to number of levels
• Nominal variables: 0 if equal, 1 otherwise

(data <- cbind(data_int, data_ord, data_nm2, data_nm3))

  int ord nm2 nm3
1   1   a   A   C
2   3   b   B   D
3   6   c   B   E

• Statistical matching with \(8 \times 10^6\) recipients vs. \(5 \times 10^4\) donors
  • Naïve approach: 20 hours
    • Compute all-pairs distances
  • Transform the problem: 20 minutes
    • Gower's -> Manhattan
    • 4 -> 6 columns
    • Use RANN.L1 for finding the \(k\) nearest donors for each recipient

devtools::install_github("krlmlr/mangow")

http://krlmlr-user15.github.io/