Match by group in R without replacement

user_9302386 0
Viewed 567

I have a data set of ~20,000 cases and ~50,000 controls. Each case has been preliminarily matched with 3 controls using a join process in SQL that asked for an exact match on some covariates (e.g., a binary variable for whether the person was enrolled in the study in 2019) and allowed a match within a range for other covariates (e.g., total days enrolled within +/- 25 days). A control could thus be listed as a possible control for more than one case.

In R, I would like to select the control for each case from each case's 3 possible matches while ensuring two things: 1) controls are not repeated, and 2) for cases for which there is no exact match, that the control is matched as closely as possible on the covariates.

I think the approach I need is to group the data by the case_id and select the control_id from each group where the difference between the covariates is minimized. The wrinkle is that I need to ensure that once a control has been selected for a case in one group, it cannot be used as the control for a second case, no matter how good the covariate match is.

Here is some sample data:

case_id <- c(1,1,1,2,2,2,3,3,3,4,4,4)
control_id <- c(5,14,7,8,9,10,11,12,13,14,15,5)
case_enrolled_2019 <- c(1,1,1,0,0,0,0,0,0,1,1,1)
control_enrolled_2019 <- c(1,1,1,0,0,0,0,0,0,1,1,1)
case_enrollment_count <- c(1,1,1,3,3,3,1,1,1,1,1,1)
control_enrollment_count <- c(2,1,1,3,3,3,1,1,1,1,1,2)
case_enrolled_days <- c(300,300,300,200,200,200,600,600,600,300,300,300)
control_enrolled_days <- c(300,300,280,200,200,210,601,610,598,300,301,300)

cbind(case_id, control_id, case_enrolled_2019, control_enrolled_2019, case_enrollment_count, control_enrollment_count, case_enrolled_days, control_enrolled_days)

I need output like:

case_id <- c(1,2,3,4)
control_id <- c(14,8,11,15)
case_enrolled_2019 <- c(1,0,0,1)
control_enrolled_2019 <- c(1,0,0,1)
case_enrollment_count <- c(1,3,1,1)
control_enrollment_count <- c(1,3,1,1)
case_enrolled_days <- c(300,200,600,300)
control_enrolled_days <- c(300,200,601,301)

cbind(case_id, control_id, case_enrolled_2019, control_enrolled_2019, case_enrollment_count, control_enrollment_count, case_enrolled_days, control_enrolled_days)
r
1 Answers

Here's a solution. It's slow because it has to loop through each case, compute the distance between it and each control unit in its block that hasn't already been matched, and then choose the nearest. Using a replicated version of this dataset with 72000 rows, it takes my computer 183 seconds. You could use parallelization to speed it up.

The first step is to create a dataset with an entry for each unique row (including control units that were initially matched with more than one treated each time) as if they were their own unit. There is probably a more elegant and generalizable way to do this, but honestly it works fine here.

 df2 <- setNames(as.data.frame(rbind(cbind(1, unique(df[,c(1,1,3,5,7)])), 
                                     cbind(0, df[,c(1,2,4,6,8)]))),
                 c("cc", "block", "id", "enrolled_2019", "enrollment_count", "enrolled_days"))

head(df2)
#>   cc block id enrolled_2019 enrollment_count enrolled_days
#> 1  1     1  1             1                1           300
#> 2  1     2  2             0                3           200
#> 3  1     3  3             0                1           600
#> 4  1     4  4             1                1           300
#> 5  0     1  5             1                2           300
#> 6  0     1 14             1                1           300

This leaves us with a dataset with 6 columns: the first cc is whether the unit is a case (1) or control (0). The second block is which initial matched block each unit is a part of. The third id is the unit's ID. The others are the variable values.

Next, we initialize the IDs of the cases and controls to match, and add a column to df2 containing whether each unit was matched or not. This is important for ensuring two rows that correspond to the same control ID are not both matched.

cases <- unique(df2$id[df2$cc == 1])
controls <- rep(NA, length(cases))
df2$matched <- FALSE

Next we loop through each case and find its match. In the loop, we first ensure there are any available control units for the case. If not, we skip the case and it remains unmatched. If so, we compute the distance between it and the remaining control units in the block using match_on from the optmatch package. This just creates a distance matrix with one row for each case and a column for each control. In this case, we have one case and at most three controls. We take the value with the minimum distance and find the control ID it corresponds to and store it in the controls vector. We also make sure to mark as "matched" any rows in df2 with the same ID so the control unit isn't reused.

for (i in seq_along(cases)) {
    if (any(df2$block==cases[i] & !df2$matched)) {
        dist <- optmatch::match_on(cc ~ enrolled_2019 + enrollment_count + enrolled_days, 
                                   data = df2[df2$block==cases[i] & !df2$matched,])
        controls[i] <- df2$id[as.numeric(colnames(dist)[which.min(dist)])]
        df2$matched[df2$id == controls[i]] <- TRUE
    }
}

Finally we extract the pairs from the matched controls. We use na.omit in case there are unmatched cases, which are discarded.

(pairs <- na.omit(data.frame(cases, controls)))
#>   cases controls
#> 1     1       14
#> 2     2        8
#> 3     3       11
#> 4     4       15

matched_df <- df[interaction(df[,1], df[,2]) %in% interaction(pairs[,1], pairs[,2]),]

Created on 2020-05-28 by the reprex package (v0.3.0)

user_63485510
Related