Replace `dplyr` grouping with fast data.table implementation, with optional strand grouping

[Arshammik]

Description

We need to replace the dplyr-based code for grouping and mutation with a faster data.table version.
The goal is to replicate the behavior of the following dplyr pipeline:

temp_eventdata_grouped <- temp_eventdata %>% 
  group_by(start_cor_id) %>% 
  mutate(start_cor_group_id = cur_group_id()) %>% 
  mutate(start_cor_group_count = n()) %>% 
  ungroup() %>% 
  group_by(end_cor_id) %>% 
  mutate(end_cor_group_id = cur_group_id()) %>% 
  mutate(end_cor_group_count = n()) %>%
  as.data.table()

Additional functionality:

• Add optional support for grouping by DNA strand (strand column) alongside start_cor_id and end_cor_id when specified by the user.
• Ensure that if the strand information is used, the groupings respect both ID and strand combinations.

Requirements

• Rewrite the above pipeline entirely using data.table for better performance.
• Keep the output structure compatible with downstream code (i.e., a data.table with the same new columns).
• Add an argument like use_strand = TRUE/FALSE to toggle whether the DNA strand should be included in grouping keys.
• Preserve the same column names:
  • start_cor_group_id
  • start_cor_group_count
  • end_cor_group_id
  • end_cor_group_count
• Avoid unnecessary copying of the data.
• Ensure the behavior when strand column is missing is either:
  • Graceful fallback (ignore strand), or
  • Explicit error (clear messaging).

Deliverables

• group_eventdata_datatable.R (main function or utility script)
• Updated documentation/comments explaining the difference between strand-aware and strand-agnostic grouping
• Unit tests to check correctness of both use_strand = TRUE and use_strand = FALSE modes
• Benchmarks comparing runtime between dplyr and data.table versions (optional)

Notes

• When use_strand = TRUE, grouping keys will be:
  • start_cor_id + strand for start_cor_group_id
  • end_cor_id + strand for end_cor_group_id
• When use_strand = FALSE, grouping will be only by start_cor_id or end_cor_id.
• data.table usage should leverage .GRP and .N efficiently:
  • .GRP gives the current group number
  • .N gives the size of the current group

Priorities

• High: Correctness of grouping logic, especially with optional strand handling.
• High: Full replacement of dplyr with data.table to improve speed.
• Medium: Defensive coding for missing or incorrectly typed strand column.
• Low: Benchmarks and performance comparison.

[Arshammik]

We replaced the original dplyr-based grouping code with a more efficient data.table implementation for assigning group IDs and counts based on coordinate identifiers:

# Use data.table for efficient grouping
if (!is.data.table(temp_eventdata)) {
  setDT(temp_eventdata)  # Convert to data.table if not already
}

# Assign group ID and count based on `start_cor_id`
temp_eventdata[, `:=`(
  start_cor_group_id = .GRP,      # Group ID for each unique start_cor_id
  start_cor_group_count = .N      # Count of entries in each group
), by = start_cor_id]

# Assign group ID and count based on `end_cor_id`
temp_eventdata[, `:=`(
  end_cor_group_id = .GRP,        # Group ID for each unique end_cor_id
  end_cor_group_count = .N        # Count of entries in each group
), by = end_cor_id]

We then validated the output by comparing the new and original versions of the multigedi_make_m1() function:

m1_new_obj <- multigedi_make_m1_new(junction_ab_object = SJ_object)
m1_old_obj <- multigedi_make_m1(junction_ab_object = SJ_object)

identical(m1_new_obj$m1_inclusion_matrix, m1_old_obj$m1_inclusion_matrix)
identical(m1_new_obj$event_data, m1_old_obj$event_data)

Output:

[1] TRUE
[1] FALSE

The inclusion matrices (m1_inclusion_matrix) are identical, confirming that downstream numeric computation is unaffected.

However, the second FALSE arises from non-identical group_id values in event_data. To isolate the issue, we excluded group_id from the comparison:

# Compare all columns except group_id
identical(
  m1_new_obj$event_data[, !("group_id"), with = FALSE],
  m1_old_obj$event_data[, !("group_id"), with = FALSE]
)

Result:

[1] TRUE

We further confirmed that the number of unique groups remains consistent between both implementations:

length(table(m1_old_obj$event_data$group_id))
length(table(m1_new_obj$event_data$group_id))

Result:

[1] 72
[1] 72

The replacement using data.table produces functionally equivalent results. The only difference is in the internal group ID labels, which are implementation-dependent but do not affect downstream computation.