fix(turnover): Fix DIANN 1.8 processing of protein turnover labeling (#140)

Author: tonywu1999

R/clean_DIANN.R

@@ -231,12 +231,17 @@
 #' \strong{not} used to validate or filter \code{ModifiedSequence}.
 #'
 #' \strong{ModifiedSequence-parsing path} (\code{has_channel = FALSE}):
-#' \code{PeptideSequence} (the retained \code{ModifiedSequence}) is matched
-#' against SILAC suffixes of the form \code{(SILAC-<AA>-H)} or
-#' \code{(SILAC-<AA>-L)}, where \code{<AA>} is any code in
-#' \code{labeledAminoAcids}.  Sequences carrying neither suffix receive
-#' \code{IsotopeLabelType = NA}.  The SILAC suffix is stripped from
-#' \code{PeptideSequence} after classification.
+#' \code{PeptideSequence} (the retained \code{ModifiedSequence}) is scanned
+#' for isotope-labeled amino acids, which appear in parentheses immediately
+#' after the labeled residue, in the form \code{(<label>-<aminoAcid>-H)} for
+#' the heavy label or \code{(<label>-<aminoAcid>-L)} for the light label.
+#' For example, \code{K(label-K-H)} marks a heavy-labeled lysine (\code{K}).
+#' Here \code{<aminoAcid>} is one of the single-letter codes in
+#' \code{labeledAminoAcids}, and \code{<label>} is the label name (e.g.
+#' \code{SILAC} or \code{label}). Sequences with no such parenthetical
+#' tags are assigned \code{IsotopeLabelType = NA}. Once classified, the
+#' parenthetical annotation is stripped out of \code{PeptideSequence},
+#' leaving the plain amino acid sequence.
 #'
 #' @param dn_input \code{data.table} after column renaming.
 #' @param labeledAminoAcids Character vector of single-letter amino acid codes
@@ -265,9 +270,9 @@
         dn_input[, Channel := NULL]
     } else {
         aa_pattern <- paste0(labeledAminoAcids, collapse = "|")
-        heavy_regex <- paste0("\\(SILAC-(?:", aa_pattern, ")-H\\)")
-        light_regex <- paste0("\\(SILAC-(?:", aa_pattern, ")-L\\)")
-        strip_regex <- paste0("\\(SILAC-(?:", aa_pattern, ")-[HL]\\)")
+        heavy_regex <- paste0("\\([^-]+-(?:", aa_pattern, ")-H\\)")
+        light_regex <- paste0("\\([^-]+-(?:", aa_pattern, ")-L\\)")
+        strip_regex <- paste0("\\([^-]+-(?:", aa_pattern, ")-[HL]\\)")
 
         dn_input <- .classifyIsotopeLabelType(dn_input, heavy_regex, light_regex)
         dn_input[, PeptideSequence := gsub(strip_regex, "", PeptideSequence, perl = TRUE)]
@@ -276,7 +281,9 @@
     if (all(is.na(dn_input[["IsotopeLabelType"]]))) {
         warning("labeledAminoAcids was provided but no rows were classified as H or L. ",
                 "Check that the input contains either a Channel column with H/L values ",
-                "or ModifiedSequence entries with (SILAC-<AA>-H)/(SILAC-<AA>-L) suffixes.")
+                "or ModifiedSequence entries with isotope-labeled amino acids", 
+                "which appear in parentheses immediately after the labeled residue", 
+                "e.g. (SILAC-K-H).")
     }
 
     dn_input

R/clean_ProteinProspector.R

@@ -2,7 +2,6 @@
 #' @param msstats_object object that inherits from MSstatsInputFiles class.
 #' @return data.table
 #' @keywords internal
-#' @noRd
 .cleanRawProteinProspector = function(msstats_object) {
     PSM = PeptideSequence = PrecursorCharge = NULL
     protein_prospector_input = getInputFile(msstats_object, "input")

R/converters_DIANNtoMSstatsFormat.R

@@ -35,12 +35,18 @@
 #'
 #' \strong{ModifiedSequence-parsing path} (DIA-NN 1.x exports without a
 #' \code{Channel} column): when \code{labeledAminoAcids} is non-\code{NULL}
-#' and no \code{Channel} column is present, each \code{ModifiedSequence} is
-#' inspected for SILAC suffixes of the form \code{(SILAC-<AA>-H)} or
-#' \code{(SILAC-<AA>-L)}, where \code{<AA>} is one of the supplied amino acid
-#' codes.  Matching sequences are classified as \code{"H"} or \code{"L"};
-#' sequences carrying neither suffix receive \code{IsotopeLabelType = NA}.
-#' The SILAC suffix is then stripped from \code{PeptideSequence}.
+#' and no \code{Channel} column is present, each \code{ModifiedSequence} is 
+#' scanned for isotope-labeled amino acids, which appear in parentheses 
+#' immediately after the labeled residue, 
+#' in the form \code{(<label>-<aminoAcid>-H)} for
+#' the heavy label or \code{(<label>-<aminoAcid>-L)} for the light label.
+#' For example, \code{K(label-K-H)} marks a heavy-labeled lysine (\code{K}).
+#' Here \code{<aminoAcid>} is one of the single-letter codes in
+#' \code{labeledAminoAcids}, and \code{<label>} is the label name (e.g.
+#' \code{SILAC} or \code{label}). Sequences with no such parenthetical
+#' tags are assigned \code{IsotopeLabelType = NA}. Once classified, the
+#' parenthetical annotation is stripped out of \code{PeptideSequence},
+#' leaving the plain amino acid sequence.
 #'
 #' When \code{NULL} (default), protein-turnover mode is disabled and all
 #' peptides receive \code{IsotopeLabelType = "Light"}.

inst/tinytest/test_clean_DIANN.R

@@ -56,7 +56,7 @@ output <- MSstatsConvert:::.cleanRawDIANN(input, MBR = FALSE, pg_qvalue_cutoff =
 expect_qvalue_cutoff(output, "GlobalPGQValue", 0.0002)
 
 # Test .assignDIANNIsotopeLabelType ---------------------------
-# Channel path: Channel column is mapped to IsotopeLabelType, then dropped
+
 dt_channel = data.table::data.table(
     PeptideSequence = c("PEPTIDEK", "PEPTIDEK", "PEPTIDEK"),
     Channel = c("H", "L", "other")
@@ -67,7 +67,6 @@ result_channel = MSstatsConvert:::.assignDIANNIsotopeLabelType(
 expect_equal(result_channel$IsotopeLabelType, c("H", "L", NA_character_))
 expect_false("Channel" %in% colnames(result_channel))
 
-# NULL path: labeledAminoAcids = NULL returns input unchanged (no IsotopeLabelType added)
 dt_null = data.table::data.table(
     PeptideSequence = c("PEPTIDEK(SILAC-K-H)", "PEPTIDEK(SILAC-K-L)")
 )
@@ -76,3 +75,27 @@ result_null = MSstatsConvert:::.assignDIANNIsotopeLabelType(
 )
 expect_false("IsotopeLabelType" %in% colnames(result_null))
 expect_equal(nrow(result_null), 2L)
+
+dt_label = data.table::data.table(
+    PeptideSequence = c("PEPTIDEK(SILAC-K-H)", "PEPTIDEK(SILAC-K-L)",
+                        "PEPTIDEK(label-K-H)", "PEPTIDEK(label-K-L)",
+                        "PEPTIDEK")
+)
+result_label = MSstatsConvert:::.assignDIANNIsotopeLabelType(
+    dt_label, labeledAminoAcids = c("K"), has_channel = FALSE
+)
+expect_equal(result_label$IsotopeLabelType, c("H", "L", "H", "L", NA_character_))
+expect_equal(unique(result_label$PeptideSequence), "PEPTIDEK")
+
+dt_multi_aa = data.table::data.table(
+    PeptideSequence = c("PEPTIDEK(SILAC-K-H)", "PEPTIDER(SILAC-R-H)",
+                        "PEPTIDEK(SILAC-K-L)", "PEPTIDER(SILAC-R-L)",
+                        "PEPTIDEAC")
+)
+result_multi_aa = MSstatsConvert:::.assignDIANNIsotopeLabelType(
+    dt_multi_aa, labeledAminoAcids = c("K", "R"), has_channel = FALSE
+)
+expect_equal(result_multi_aa$IsotopeLabelType,
+             c("H", "H", "L", "L", NA_character_))
+expect_equal(sort(unique(result_multi_aa$PeptideSequence)),
+             c("PEPTIDEAC", "PEPTIDEK", "PEPTIDER"))