FlexRL-RL-block-vignette
FlexRL is a package proposing flexible probabilistic
Record Linkage with diagnostic tools
for downstream inference on linked data. The record linkage algorithm
StEM() uses a Stochastic
Expectation Maximisation approach to
combine 2 data sources and outputs the set of records referring to the
same entities.
More details on the record linkage methodology; on the false discovery proportion estimation; and on the potential risks of downstream inference on linked data (Robach et al., in preparation).
Several other algorithms for record linkage are available in R. The
link_with_*() wrappers of FlexRL run them through one
common interface, so that the same diagnostics and false discovery
proportion (FDP) estimators apply to all of them:
link_with_FlexRL(), link_with_BRL(),
link_with_fastLink(), link_with_multilink(),
link_with_reclin2(), link_with_diyar(),
link_with_fedmatch(). Each wrapper takes the two encoded
sources and a named list of arguments forwarded to the underlying
package, and returns a list with idxA, idxB
(row indices of the linked records in dataA and dataB) and
LinkScore (the linkage scores, NULL when the
package does not return any).
This vignette runs these wrappers on the same pipeline for: - a regional subset of the NLTCS data (1982 and 1994 waves) - a regional subset of the SHIW data (2016 and 2020 waves) - the full SHIW data with and without regional blocking
One pipeline for all the wrappers
wrapper_args <- function(method, prep_data, ...) {
PIVs <- names(prep_data$PIVs_config)
args <- c(...)
switch(method,
FlexRL = c(list(data = prep_data),
args[intersect(names(args), names(formals(FlexRL::StEM)))]),
BRL = c(list(flds = PIVs, types = rep("bi", length(PIVs)), nIter = 1000),
args[intersect(names(args), c(names(formals(BRL::compareRecords)), names(formals(BRL::bipartiteGibbs))))]),
fastLink = c(list(varnames = PIVs, tol.em = 1e-06, threshold.match = 0.5, return.all = TRUE, n.cores = 1),
args[intersect(names(args), names(formals(fastLink::fastLink)))]),
multilink = c(list(types = rep("bi", length(PIVs)), breaks = stats::setNames(as.list(rep(NA, length(PIVs))), PIVs),
duplicates = c(0, 0), n_iter = 200, max_cc_size = 50),
args[intersect(names(args), c(names(formals(multilink::create_comparison_data)),
names(formals(multilink::reduce_comparison_data)),
names(formals(multilink::specify_prior)),
names(formals(multilink::gibbs_sampler)),
names(formals(multilink::find_bayes_estimate)),
names(formals(multilink::relabel_bayes_estimate))))]),
reclin2 = c(list(on = PIVs, formula = stats::as.formula(paste("~", paste(PIVs, collapse = "+"))),
type = "mpost", add = TRUE, variable = "selected", score = "mpost", threshold = 0.5),
args[intersect(names(args), c(names(formals(reclin2::pair)),
names(formals(reclin2::compare_pairs)),
names(formals(reclin2::problink_em)),
names(formals(stats::predict)),
names(formals(reclin2::select_threshold))))]),
diyar = c(list(probabilistic = TRUE, return_weights = TRUE),
args[intersect(names(args), c(names(formals(diyar::prob_score_range)),
names(formals(diyar::links_wf_probabilistic))))]),
fedmatch = c(list(by = PIVs, match_type = "multivar", unique_key_1 = "unique_key_A",
unique_key_2 = "unique_key_B", multivar_settings = fedmatch::build_multivar_settings(
compare_type = rep("indicator", length(PIVs)),
wgts = rep(1 / length(PIVs), length(PIVs)))),
args[intersect(names(args), names(formals(fedmatch::merge_plus)))])
)
}
run_wrapper <- function(method, prep_data, arguments) {
PIVs <- names(prep_data$PIVs_config)
link_with <- get(paste0("link_with_", method), envir = asNamespace("FlexRL"))
if (method == "diyar") {
prep_data$encodedA <- prep_data$encodedA[, c(PIVs, "local_id", "source")]
prep_data$encodedB <- prep_data$encodedB[, c(PIVs, "local_id", "source")]
}
tryCatch(link_with(prep_data$encodedA, prep_data$encodedB, arguments),
error = function(e) {
if (grepl("cannot allocate|memory|bad_alloc|too large|long vectors", conditionMessage(e), ignore.case = TRUE)) {
message(method, ": memory error (", conditionMessage(e), "), results reported as NA")
return(NULL)
}
stop(e)
})
}
evaluate_linkage <- function(fit, true_pairs, threshold = 0.5) {
if (is.null(fit)) return(rep(NA_real_, 5))
delta_result <- data.frame(idxA = fit$idxA, idxB = fit$idxB)
if (!is.null(fit$LinkScore)) delta_result <- delta_result[fit$LinkScore > threshold, , drop = FALSE]
linked_pairs <- do.call(paste, c(delta_result, list(sep = "_")))
true_positive <- length(intersect(linked_pairs, true_pairs))
false_positive <- length(setdiff(linked_pairs, true_pairs))
false_negative <- length(setdiff(true_pairs, linked_pairs))
c(TP = true_positive, FP = false_positive, FN = false_negative,
sensitivity = true_positive / (true_positive + false_negative),
FDP = false_positive / (true_positive + false_positive))
}
results_rows <- c("TP",
"FP",
"FN",
"sensitivity",
"FDP",
"hat FDP score",
"hat FDP score*",
"hat FDP synth*",
"max SMD",
"min support IoU",
"max MMD")
run_all_wrappers <- function(prep_data, PIVs_type, methods, FDP_estimation = TRUE,
maxIter4CV = 3, n_repeats = 5, ...) {
args <- list(...)
op <- options()
on.exit(options(op))
PIVs <- names(prep_data$PIVs_config)
true_pairs <- do.call(paste, c(prep_data$true_pairs, list(sep = "_")))
df_results <- data.frame(matrix(NA, nrow = length(results_rows), ncol = 0))
rownames(df_results) <- results_rows
naive_fit <- naive_linkage(PIVs = PIVs, encodedA = prep_data$encodedA, encodedB = prep_data$encodedB)
df_results[1:5, "Naive"] <- evaluate_linkage(list(idxA = naive_fit[[1]], idxB = naive_fit[[2]]), true_pairs)
for (method in methods) {
arguments <- wrapper_args(method, prep_data, args)
fit <- run_wrapper(method, prep_data, arguments)
df_results[1:5, method] <- evaluate_linkage(fit, true_pairs)
if (is.null(fit) || length(fit$idxA) == 0) next
encodedA <- prep_data$encodedA
encodedB <- prep_data$encodedB
if (method == "diyar") {
encodedA <- encodedA[, c(PIVs, "local_id", "source")]
encodedB <- encodedB[, c(PIVs, "local_id", "source")]
}
diag <- do.call(RL_diagnostics, c(list(fit, encodedA, encodedB, PIVs, PIVs_type, PIVs,
true_pairs = prep_data$true_pairs, FDP_estimation = FDP_estimation,
RL_method = method, data = prep_data, maxIter4CV = maxIter4CV, n_repeats = n_repeats),
arguments))
if (FDP_estimation) {
df_results[6, method] <- diag$FDP_measures$curve$RL_FDP_score[1]
df_results[7, method] <- diag$FDP_measures$curve$augmRL_FDP_score[1]
df_results[8, method] <- diag$FDP_measures$curve$augmRL_FDP_synth[1]
}
df_results[9, method] <- max(diag$discrepancy_measures[1, grep("smd", names(diag$discrepancy_measures))])
df_results[10, method] <- min(diag$discrepancy_measures[1, grep("iou", names(diag$discrepancy_measures))])
df_results[11, method] <- max(diag$discrepancy_measures[1, grep("mmd", names(diag$discrepancy_measures))])
}
df_results
}
show_results <- function(df_results) {
cat("Record Linkage\n"); print(df_results[1:5, ], 2)
cat("\nFDP estimation\n"); print(df_results[6:8, ], 2)
cat("\nLinked divergence\n"); print(df_results[9:11, ], 2)
}
methods <- c("FlexRL", "BRL", "multilink", "fastLink", "reclin2", "diyar", "fedmatch")We will fill a table of results throughout the vignette, gathering for each record linkage method: linkage performance, data quality diagnostics: FDP estimation, measures of divergence of linked data from sources.
Regional subset of the NLTCS data
We use 4 Partially Identifying
Variables to link the records: sex,
birth month and birth year are stable, the state code may
change over time but we cannot model these changes so we set this
variable as flexible.
df1982 <- read.csv("ex-reg-NLTCS-82.csv", row.names = 1)
df1994 <- read.csv("ex-reg-NLTCS-94.csv", row.names = 1)
PIVs_config <- list(
SEX = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)),
DOB_MONTH = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)),
DOB_YEAR = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)),
STATE = list(dynamics = "flexible", bound_mistakes = c(NA, NA), fix_mistakes = c(NA, NA)))
PIVs <- names(PIVs_config)
PIVs_type <- list(SEX = FALSE, DOB_MONTH = TRUE, DOB_YEAR = TRUE, STATE = FALSE)prepare_data() removes the records for which the linking
variables are outside of the common support, checks that PIVs are not
strongly associated with each other (which may degrade record linkage
performance), maps PIV values to natural numbers and encodes missing
values as 0. The StEM algorithm denotes the smallest source as
A and the largest as B. The true pairs are
known from the identifier ID (they are indexed on the rows of
encodedA and encodedB, after the support
restriction).
prep_data <- prepare_data(df1982, df1994, "1982", "1994", PIVs_config, same_mistakes = TRUE, uniq_id = "ID")
#> '1982' is the larger source, saved as B; '1994' saved as A.
c(n_A = nrow(prep_data$encodedA), n_B = nrow(prep_data$encodedB), n_links = nrow(prep_data$true_pairs))
#> n_A n_B n_links
#> 491 901 332
prep_data$n_values
#> SEX DOB_MONTH DOB_YEAR STATE
#> 2 12 28 9show_results(df_results_NLTCS)
#> Record Linkage
#> Naive FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> TP 320.00 221.00 212.000 236.00 255.00 320.00 126.000 266.00
#> FP 163.00 28.00 15.000 44.00 88.00 163.00 5.000 223.00
#> FN 12.00 111.00 120.000 96.00 77.00 12.00 206.000 66.00
#> sensitivity 0.96 0.67 0.639 0.71 0.77 0.96 0.380 0.80
#> FDP 0.34 0.11 0.066 0.16 0.26 0.34 0.038 0.46
#>
#> FDP estimation
#> Naive FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> hat FDP score NA 0.16 0.151 NaN 0.30 0.49 NaN 0.075
#> hat FDP score* NA 0.17 0.169 NaN 0.31 NA NaN 0.071
#> hat FDP synth* NA 0.13 0.093 0.17 0.22 NA 0.12 NaN
#>
#> Linked divergence
#> Naive FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> max SMD NA 0.210 0.149 0.141 0.141 0.177 0.43 0.1125
#> min support IoU NA 0.875 0.875 0.875 0.875 0.875 0.88 0.8750
#> max MMD NA 0.015 0.014 0.014 0.016 0.025 0.03 0.0093Regional subset of the SHIW data
We use 4 Partially Identifying
Variables to link the records: sex and
birth year are stable, education level and marital status
are flexible.
df2016 <- read.csv("ex-reg-SHIW-16.csv", row.names = 1)
df2020 <- read.csv("ex-reg-SHIW-20.csv", row.names = 1)
PIVs_config <- list(
ANASCI = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)),
SESSO = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)),
STUDIO = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)),
STACIV = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA)))
PIVs <- names(PIVs_config)
PIVs_type <- list(ANASCI = TRUE, SESSO = FALSE, STUDIO = FALSE, STACIV = FALSE)
prep_data <- prepare_data(df2016, df2020, "2016", "2020", PIVs_config, same_mistakes = TRUE, uniq_id = "ID")
#> '2016' is the larger source, saved as B; '2020' saved as A.
c(n_A = nrow(prep_data$encodedA), n_B = nrow(prep_data$encodedB), n_links = nrow(prep_data$true_pairs))
#> n_A n_B n_links
#> 42 49 27
prep_data$n_values
#> ANASCI SESSO STUDIO STACIV
#> 28 2 5 4show_results(df_results_SHIW)
#> Record Linkage
#> Naive FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> TP 20.00 17.00 21.00 23.00 19.00 27.00 13.000 19.00
#> FP 8.00 7.00 14.00 14.00 11.00 52.00 1.000 20.00
#> FN 7.00 10.00 6.00 4.00 8.00 0.00 14.000 8.00
#> sensitivity 0.74 0.63 0.78 0.85 0.70 1.00 0.481 0.70
#> FDP 0.29 0.29 0.40 0.38 0.37 0.66 0.071 0.51
#>
#> FDP estimation
#> Naive FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> hat FDP score NA 0.14 0.20 NaN 0.10 0.14 NaN 0.09
#> hat FDP score* NA 0.14 0.21 NaN 0.12 0.14 NaN 0.09
#> hat FDP synth* NA 0.10 0.21 0.41 0.16 0.47 0.19 NaN
#>
#> Linked divergence
#> Naive FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> max SMD NA 0.357 0.081 0.120 0.190 0.123 0.64 0.168
#> min support IoU NA 0.333 1.000 1.000 1.000 1.000 0.83 1.000
#> max MMD NA 0.037 0.019 0.018 0.028 0.018 0.14 0.027Full SHIW data
df2016 <- read.csv("SHIW2016.csv", row.names = 1)
df2020 <- read.csv("SHIW2020.csv", row.names = 1)
df2016 <- df2016[!df2016$ID %in% df2016$ID[duplicated(df2016$ID)], ]
df2020 <- df2020[!df2020$ID %in% df2020$ID[duplicated(df2020$ID)], ]
PIVs_config_full <- c(PIVs_config,
list(IREG = list(dynamics = "stable", bound_mistakes = c(0.10, 0.10), fix_mistakes = c(NA, NA))))
PIVs_full <- names(PIVs_config_full)
PIVs_type_full <- c(PIVs_type, IREG = FALSE)
prep_data_full <- prepare_data(df2016, df2020, "2016", "2020", PIVs_config_full, same_mistakes = TRUE, uniq_id = "ID")
#> '2016' is the larger source, saved as B; '2020' saved as A.
true_pairs_full <- do.call(paste, c(prep_data_full$true_pairs, list(sep = "_")))df_results_full <- data.frame(matrix(NA, nrow = 6, ncol = 0))
rownames(df_results_full) <- c("TP", "FP", "FN", "sensitivity", "FDP", "minutes")
for (method in methods) {
t0 <- Sys.time()
args <- list(StEM_iter = 20, StEM_burnin = 10, gibbs_iter = 20, gibbs_burnin = 10,
n_post_sample = 100, n_iter = 500, max_cc_size = 50, n.cores = 4)
fit <- run_wrapper(method, prep_data_full, wrapper_args(method, prep_data_full, args))
df_results_full[1:5, method] <- evaluate_linkage(fit, true_pairs_full)
df_results_full[6, method] <- if (is.null(fit)) NA else round(as.numeric(Sys.time() - t0, units = "mins"))
}
# on Apple M1 Max ARM chip, 8-core CPU, 32 GB of memory
# fastLink: memory error / multilink: memory error / BRL: memory error / diyar: memory error
# we reduce the number of iterations of FlexRL to gain in computational timecat("Record Linkage performance on large files")
#> Record Linkage performance on large files
round(df_results_full,2)
#> FlexRL BRL multilink fastLink reclin2 diyar fedmatch
#> TP 1046.00 NA NA NA 5936.00 NA 2756.00
#> FP 1560.00 NA NA NA 1917562.00 NA 12161.00
#> FN 5384.00 NA NA NA 494.00 NA 3674.00
#> sensitivity 0.16 NA NA NA 0.92 NA 0.43
#> FDP 0.60 NA NA NA 1.00 NA 0.82
#> minutes 107.00 NA NA NA 2.00 NA 7.00Note that when methods do not ensure the one-to-one assignment constraint the number of False Positives may be high due to many-to-many falsely linked pairs of records.
Loop over the regions of the SHIW data with blocking
Blocking on the region code IREG splits the task into 20
smaller ones. Within a region, IREG is constant and the
linkage uses the 4 other PIVs. Region 2 is the regional subset used
above.
regions <- sort(unique(c(df2016$IREG, df2020$IREG)))
results_IREG <- list()
for (r in regions) {
block2016 <- df2016[df2016$IREG == r, ]
block2020 <- df2020[df2020$IREG == r, ]
prep_data_r <- prepare_data(block2016, block2020, "2016", "2020", PIVs_config, same_mistakes = TRUE, uniq_id = "ID")
true_pairs_r <- do.call(paste, c(prep_data_r$true_pairs, list(sep = "_")))
for (method in methods) {
t0 <- Sys.time()
fit <- run_wrapper(method, prep_data_r, wrapper_args(method, prep_data_r))
results_IREG[[length(results_IREG) + 1]] <- data.frame(
IREG = r, n_A = nrow(prep_data_r$encodedA),
n_B = nrow(prep_data_r$encodedB), n_links = length(true_pairs_r),
method = method, round(t(evaluate_linkage(fit, true_pairs_r)),2),
seconds = round(as.numeric(Sys.time() - t0, units = "secs")))
}
}
results_IREG <- do.call(rbind, results_IREG)cat("Record Linkage on regional subsets")
#> Record Linkage on regional subsets
results_IREG
#> IREG n_A n_B n_links method TP FP FN sensitivity FDP seconds
#> 1 1 965 1272 388 FlexRL 35 46 353 0.09 0.57 1181
#> 2 1 965 1272 388 BRL 21 13 367 0.05 0.38 13
#> 3 1 965 1272 388 multilink 50 243 338 0.13 0.83 168
#> 4 1 965 1272 388 fastLink 133 654 255 0.34 0.83 11
#> 5 1 965 1272 388 reclin2 350 9033 38 0.90 0.96 0
#> 6 1 965 1272 388 diyar 36 38 352 0.09 0.51 4
#> 7 1 965 1272 388 fedmatch 135 830 253 0.35 0.86 9
#> 8 2 42 49 27 FlexRL 16 10 11 0.59 0.38 17
#> 9 2 42 49 27 BRL 19 14 8 0.70 0.42 0
#> 10 2 42 49 27 multilink 23 14 4 0.85 0.38 0
#> 11 2 42 49 27 fastLink 19 11 8 0.70 0.37 1
#> 12 2 42 49 27 reclin2 27 52 0 1.00 0.66 0
#> 13 2 42 49 27 diyar 13 1 14 0.48 0.07 0
#> 14 2 42 49 27 fedmatch 19 20 8 0.70 0.51 0
#> 15 3 1722 1791 636 FlexRL 42 33 594 0.07 0.44 1855
#> 16 3 1722 1791 636 BRL 13 1 623 0.02 0.07 31
#> 17 3 1722 1791 636 multilink 37 240 599 0.06 0.87 625
#> 18 3 1722 1791 636 fastLink 220 1096 416 0.35 0.83 20
#> 19 3 1722 1791 636 reclin2 616 26566 20 0.97 0.98 1
#> 20 3 1722 1791 636 diyar 40 26 596 0.06 0.39 9
#> 21 3 1722 1791 636 fedmatch 300 1422 336 0.47 0.83 17
#> 22 4 445 693 254 FlexRL 60 48 194 0.24 0.44 483
#> 23 4 445 693 254 BRL 52 27 202 0.20 0.34 3
#> 24 4 445 693 254 multilink 71 105 183 0.28 0.60 24
#> 25 4 445 693 254 fastLink 117 250 137 0.46 0.68 4
#> 26 4 445 693 254 reclin2 247 2784 7 0.97 0.92 0
#> 27 4 445 693 254 diyar 13 7 241 0.05 0.35 2
#> 28 4 445 693 254 fedmatch 151 294 103 0.59 0.66 4
#> 29 5 1042 1122 446 FlexRL 57 43 389 0.13 0.43 1209
#> 30 5 1042 1122 446 BRL 37 9 409 0.08 0.20 10
#> 31 5 1042 1122 446 multilink 68 318 378 0.15 0.82 147
#> 32 5 1042 1122 446 fastLink 184 711 262 0.41 0.79 8
#> 33 5 1042 1122 446 reclin2 435 11856 11 0.98 0.96 0
#> 34 5 1042 1122 446 diyar 29 11 417 0.07 0.28 5
#> 35 5 1042 1122 446 fedmatch 234 807 212 0.52 0.78 10
#> 36 6 326 359 165 FlexRL 52 24 113 0.32 0.32 275
#> 37 6 326 359 165 BRL 54 29 111 0.33 0.35 1
#> 38 6 326 359 165 multilink 65 68 100 0.39 0.51 8
#> 39 6 326 359 165 fastLink 108 147 57 0.65 0.58 1
#> 40 6 326 359 165 reclin2 164 1152 1 0.99 0.88 0
#> 41 6 326 359 165 diyar 22 6 143 0.13 0.21 1
#> 42 6 326 359 165 fedmatch 86 240 79 0.52 0.74 2
#> 43 7 400 494 163 FlexRL 46 37 117 0.28 0.45 353
#> 44 7 400 494 163 BRL 42 28 121 0.26 0.40 2
#> 45 7 400 494 163 multilink 51 66 112 0.31 0.56 18
#> 46 7 400 494 163 fastLink 93 209 70 0.57 0.69 2
#> 47 7 400 494 163 reclin2 159 1625 4 0.98 0.91 0
#> 48 7 400 494 163 diyar 45 27 118 0.28 0.38 1
#> 49 7 400 494 163 fedmatch 91 304 72 0.56 0.77 3
#> 50 8 1120 1217 418 FlexRL 41 37 377 0.10 0.47 1221
#> 51 8 1120 1217 418 BRL 27 23 391 0.06 0.46 13
#> 52 8 1120 1217 418 multilink 66 233 352 0.16 0.78 209
#> 53 8 1120 1217 418 fastLink 163 617 255 0.39 0.79 9
#> 54 8 1120 1217 418 reclin2 414 12047 4 0.99 0.97 0
#> 55 8 1120 1217 418 diyar 26 18 392 0.06 0.41 5
#> 56 8 1120 1217 418 fedmatch 165 955 253 0.39 0.85 11
#> 57 9 900 1126 412 FlexRL 75 76 337 0.18 0.50 1207
#> 58 9 900 1126 412 BRL 68 37 344 0.17 0.35 10
#> 59 9 900 1126 412 multilink 118 250 294 0.29 0.68 116
#> 60 9 900 1126 412 fastLink 178 577 234 0.43 0.76 9
#> 61 9 900 1126 412 reclin2 400 8718 12 0.97 0.96 0
#> 62 9 900 1126 412 diyar 52 35 360 0.13 0.40 4
#> 63 9 900 1126 412 fedmatch 187 710 225 0.45 0.79 9
#> 64 10 592 614 316 FlexRL 92 57 224 0.29 0.38 680
#> 65 10 592 614 316 BRL 84 29 232 0.27 0.26 3
#> 66 10 592 614 316 multilink 113 141 203 0.36 0.56 34
#> 67 10 592 614 316 fastLink 166 299 150 0.53 0.64 4
#> 68 10 592 614 316 reclin2 313 3620 3 0.99 0.92 0
#> 69 10 592 614 316 diyar 60 20 256 0.19 0.25 2
#> 70 10 592 614 316 fedmatch 161 426 155 0.51 0.73 5
#> 71 11 439 686 278 FlexRL 87 50 191 0.31 0.36 575
#> 72 11 439 686 278 BRL 75 30 203 0.27 0.29 3
#> 73 11 439 686 278 multilink 104 118 174 0.37 0.53 22
#> 74 11 439 686 278 fastLink 156 224 122 0.56 0.59 4
#> 75 11 439 686 278 reclin2 277 2851 1 1.00 0.91 0
#> 76 11 439 686 278 diyar 33 15 245 0.12 0.31 2
#> 77 11 439 686 278 fedmatch 142 297 136 0.51 0.68 4
#> 78 12 959 1115 379 FlexRL 48 32 331 0.13 0.40 885
#> 79 12 959 1115 379 BRL 43 18 336 0.11 0.30 11
#> 80 12 959 1115 379 multilink 64 182 315 0.17 0.74 152
#> 81 12 959 1115 379 fastLink 163 459 216 0.43 0.74 7
#> 82 12 959 1115 379 reclin2 369 8630 10 0.97 0.96 0
#> 83 12 959 1115 379 diyar 46 25 333 0.12 0.35 4
#> 84 12 959 1115 379 fedmatch 215 744 164 0.57 0.78 9
#> 85 13 315 420 141 FlexRL 41 37 100 0.29 0.47 296
#> 86 13 315 420 141 BRL 36 27 105 0.26 0.43 1
#> 87 13 315 420 141 multilink 50 56 91 0.35 0.53 9
#> 88 13 315 420 141 fastLink 75 164 66 0.53 0.69 1
#> 89 13 315 420 141 reclin2 137 1103 4 0.97 0.89 0
#> 90 13 315 420 141 diyar 8 4 133 0.06 0.33 1
#> 91 13 315 420 141 fedmatch 80 235 61 0.57 0.75 2
#> 92 14 181 212 118 FlexRL 57 16 61 0.48 0.22 144
#> 93 14 181 212 118 BRL 59 19 59 0.50 0.24 0
#> 94 14 181 212 118 multilink 70 36 48 0.59 0.34 2
#> 95 14 181 212 118 fastLink 82 63 36 0.69 0.43 1
#> 96 14 181 212 118 reclin2 117 510 1 0.99 0.81 0
#> 97 14 181 212 118 diyar 32 4 86 0.27 0.11 1
#> 98 14 181 212 118 fedmatch 69 110 49 0.58 0.61 1
#> 99 15 1324 1734 784 FlexRL 91 76 693 0.12 0.46 2076
#> 100 15 1324 1734 784 BRL 62 32 722 0.08 0.34 22
#> 101 15 1324 1734 784 multilink 128 446 656 0.16 0.78 294
#> 102 15 1324 1734 784 fastLink 232 902 552 0.30 0.80 15
#> 103 15 1324 1734 784 reclin2 746 19909 38 0.95 0.96 0
#> 104 15 1324 1734 784 diyar 27 11 757 0.03 0.29 9
#> 105 15 1324 1734 784 fedmatch 274 1050 510 0.35 0.79 14
#> 106 16 1033 1041 465 FlexRL 61 52 404 0.13 0.46 1172
#> 107 16 1033 1041 465 BRL 38 21 427 0.08 0.36 8
#> 108 16 1033 1041 465 multilink 92 335 373 0.20 0.78 131
#> 109 16 1033 1041 465 fastLink 173 666 292 0.37 0.79 8
#> 110 16 1033 1041 465 reclin2 444 10092 21 0.95 0.96 0
#> 111 16 1033 1041 465 diyar 13 5 452 0.03 0.28 6
#> 112 16 1033 1041 465 fedmatch 233 800 232 0.50 0.77 10
#> 113 17 155 235 103 FlexRL 54 22 49 0.52 0.29 127
#> 114 17 155 235 103 BRL 57 24 46 0.55 0.30 0
#> 115 17 155 235 103 multilink 64 52 39 0.62 0.45 2
#> 116 17 155 235 103 fastLink 69 58 34 0.67 0.46 1
#> 117 17 155 235 103 reclin2 98 375 5 0.95 0.79 0
#> 118 17 155 235 103 diyar 39 5 64 0.38 0.11 0
#> 119 17 155 235 103 fedmatch 57 98 46 0.55 0.63 1
#> 120 18 324 394 174 FlexRL 53 28 121 0.30 0.34 306
#> 121 18 324 394 174 BRL 57 27 117 0.33 0.32 1
#> 122 18 324 394 174 multilink 70 68 104 0.40 0.49 9
#> 123 18 324 394 174 fastLink 97 150 77 0.56 0.61 1
#> 124 18 324 394 174 reclin2 171 1270 3 0.98 0.88 0
#> 125 18 324 394 174 diyar 47 12 127 0.27 0.20 1
#> 126 18 324 394 174 fedmatch 90 233 84 0.52 0.72 2
#> 127 19 1213 1353 414 FlexRL 36 47 378 0.09 0.57 1408
#> 128 19 1213 1353 414 BRL 15 8 399 0.04 0.35 15
#> 129 19 1213 1353 414 multilink 53 304 361 0.13 0.85 237
#> 130 19 1213 1353 414 fastLink 138 624 276 0.33 0.82 11
#> 131 19 1213 1353 414 reclin2 397 14624 17 0.96 0.97 0
#> 132 19 1213 1353 414 diyar 22 22 392 0.05 0.50 7
#> 133 19 1213 1353 414 fedmatch 160 1047 254 0.39 0.87 12
#> 134 20 709 774 332 FlexRL 87 50 245 0.26 0.36 735
#> 135 20 709 774 332 BRL 76 30 256 0.23 0.28 5
#> 136 20 709 774 332 multilink 117 178 215 0.35 0.60 60
#> 137 20 709 774 332 fastLink 165 317 167 0.50 0.66 4
#> 138 20 709 774 332 reclin2 323 4861 9 0.97 0.94 0
#> 139 20 709 774 332 diyar 22 11 310 0.07 0.33 3
#> 140 20 709 774 332 fedmatch 209 498 123 0.63 0.70 6Averaged over the regions and total computing time:
summary_IREG <- do.call(rbind, lapply(split(results_IREG, results_IREG$method), function(d) {
data.frame(method = d$method[1],
FDP = mean(d$FDP, na.rm=TRUE),
sensitivity = mean(d$sensitivity, na.rm=TRUE),
TP = mean(d$TP, na.rm=TRUE), FP = mean(d$FP, na.rm=TRUE), FN = mean(d$FN, na.rm=TRUE),
minutes = sum(d$seconds) / 60)
}))
rownames(summary_IREG) <- summary_IREG$method
round(summary_IREG[,2:7],2)
#> FDP sensitivity TP FP FN minutes
#> BRL 0.32 0.23 46.75 22.30 273.90 2.53
#> diyar 0.30 0.15 31.25 15.15 289.40 1.12
#> fastLink 0.68 0.49 136.55 409.90 184.10 2.03
#> fedmatch 0.74 0.51 152.90 556.00 167.75 2.18
#> FlexRL 0.42 0.25 56.55 41.05 264.10 270.08
#> multilink 0.63 0.32 73.70 172.65 246.95 37.78
#> reclin2 0.91 0.97 310.20 7083.90 10.45 0.02colours <- grDevices::hcl.colors(length(methods), "Dark 3")
names(colours) <- methods
plot(NA, xlim = c(0, 1), ylim = c(0, 1), xlab = "sensitivity", ylab = "FDP",
main = "SHIW linked region by region (one point per region)")
for (method in methods) {
d <- results_IREG[results_IREG$method == method, ]
points(d$sensitivity, d$FDP, col = colours[method], pch = 19)
}
legend("bottomright", legend = methods, col = colours, pch = 19, bty = "n")
plot of chunk shiw-regional-blocking-plot
The 4 PIVs have little discriminating power on these data (birth
years, sex, 6 education levels, 4 marital statuses) and the two waves
overlap on about 40% of their records (6.430 true links, 14.917 records
in A, 16.445 records in B), so the task is hard for every package.
fastLink, multilink, BRL and
diyar are not scalable (out of memory errors).
We compare methods with blocking over the 20 regions:
reclin2 (EM on the comparison patterns, threshold at 0.5 on
the posterior) links a lot of records reaching a high sensitivity at the
price of a high FDP; fastLink, fedmatch and
multilink link about half of the true pairs with a FDP
around 0.7; BRL and diyar are the most
conservative, with a FDP around 0.3 and a sensitivity around 0.2;
FlexRL performance is between these conservative and
permissive methods.
In terms of computational time however, FlexRL and
multilink are more than 10 times slower than the other
methods. The results were computed with the code above, on a machine
with Apple M1 Max ARM chip, 8-core CPU, 32 GB of memory.