Address issue #67 by adding parameter "protein_modlist_file" which points to a file that specifies restricting variable mods to apply only to a subset of proteins. Contents of the text file pointed to by that parameter are modification number and an accession string, one per line. This functionality current only works with standard Comet; support for FI will be addressed next. Also simplify WriteVariableMod() in both CometWritePepXML and CometWriteMzIdentML.

Author: jke000

Comet.cpp

@@ -366,6 +366,13 @@ void LoadParameters(char *pszParamsFile,
                sscanf(szParamVal, "%255s", szTxtFileExt);
                pSearchMgr->SetParam("text_file_extension", szTxtFileExt, szTxtFileExt);
             }
+            else if (!strcmp(szParamName, "protein_modslist_file"))
+            {
+               char szTmp[256];
+               szTmp[0] = '\0';
+               sscanf(szParamVal, "%255s", szTmp);
+               pSearchMgr->SetParam("protein_modslist_file", szTmp, szTmp);
+            }
             else if (!strcmp(szParamName, "explicit_deltacn"))
             {
                sscanf(szParamVal, "%d", &iIntParam);
@@ -1758,8 +1765,14 @@ variable_mod15 = 0.0 X 0 3 -1 0 0 0.0\n");
 
    fprintf(fp,
 "max_variable_mods_in_peptide = 5\n\
-require_variable_mod = 0\n\
-\n\
+require_variable_mod = 0\n");
+   if (iPrintParams == 2)
+   {
+      fprintf(fp, "protein_modslist_file =                # limit variable mods to subset of specified proteins if this file is specified & present\n");
+   }
+
+   fprintf(fp,
+"\n\
 #\n\
 # fragment ions\n\
 #\n\

CometSearch/CometDataInternal.h

@@ -130,7 +130,6 @@ struct Options
    int bClipNtermMet;            // 0=leave protein sequences alone; 1=also consider w/o N-term methionine
    int bClipNtermAA;             // 0=leave peptide sequences as-is; 1=clip N-term amino acid from every peptide
    int bSkipAlreadyDone;         // 0=search everything; 1=don't re-search if .out exists
-// int bSkipUpdateCheck;         // 0=do not check for updates; 1=check for updates
    int bMango;                   // 0=normal; 1=Mango x-link ms2 input
    int bScaleFragmentNL;         // 0=no; 1=scale fragment NL for each modified residue contained in fragment
    int bCreateIndex;             // 0=normal search; 1=create peptide index file
@@ -191,7 +190,6 @@ struct Options
       bClipNtermMet = a.bClipNtermMet;
       bClipNtermAA = a.bClipNtermAA;
       bSkipAlreadyDone = a.bSkipAlreadyDone;
-//    bSkipUpdateCheck = a.bSkipUpdateCheck;
       bMango = a.bMango;
       bScaleFragmentNL = a.bScaleFragmentNL;
       bCreateIndex = a.bCreateIndex;
@@ -566,6 +564,8 @@ struct VarModParams
    int     iMaxPermutations;
    VarMods varModList[VMODS];
    char    cModCode[VMODS];          // mod characters
+   string  sProteinLModsListFile;                 // file containing list of proteins to restrict application of varmods to
+   multimap<int, string> mmapProteinLModsList;    // <varmod#, protein name> vector read from sProteinModsListFile if present
 
    VarModParams& operator=(VarModParams& a)
    {
@@ -822,6 +822,14 @@ struct StaticParams
       variableModParameters.cModCode[6] = '%';
       variableModParameters.cModCode[7] = '!';
       variableModParameters.cModCode[8] = '+';
+      for (int i = 9; i < VMODS; ++i)
+      {
+         int iAscii = 88 + i;    //start with lower case 'a' ASCII 97
+         if (iAscii <= 125)      // thru '}' which is ASCII 125
+            variableModParameters.cModCode[i] = (char)(iAscii);
+         else
+            variableModParameters.cModCode[i] = '_';
+      }
 
       variableModParameters.iMaxVarModPerPeptide = 5;
       variableModParameters.iMaxPermutations = MAX_PERMUTATIONS;
@@ -866,7 +874,6 @@ struct StaticParams
       options.bResolveFullPaths = 1;
 
       options.bSkipAlreadyDone = 1;
-//    options.bSkipUpdateCheck = 0;
       options.bMango = 0;
       options.bScaleFragmentNL = 0;
       options.bCreateIndex = 0;

CometSearch/CometSearch.cpp

@@ -23,6 +23,7 @@
 #include "ModificationsPermuter.h"
 
 #include <stdio.h>
+#include <string.h>
 #include <sstream>
 #include <bitset>
 
@@ -1569,6 +1570,7 @@ bool CometSearch::SearchForPeptides(struct sDBEntry dbe,
    int iStartPos = 0;
    int iEndPos = 0;
    int piVarModCounts[VMODS];
+   bool pbVarModProteinFilter[VMODS];  // default true; set to false if a mod has a protein filter that does not match this protein
    int iWhichIonSeries;
    int ctIonSeries;
    int ctLen;
@@ -1612,6 +1614,41 @@ bool CometSearch::SearchForPeptides(struct sDBEntry dbe,
    if (g_staticParams.options.bClipNtermAA) // skip the N-term residue of every peptide
       iStartPos = 1;
 
+   for (int i = 0; i < VMODS; ++i)
+   {
+      piVarModCounts[i] = 0;
+      pbVarModProteinFilter[i] = true;
+   }
+
+   // If variable modifications protein filter is applied, check if current sequence
+   // is in the protein filter list.  Check each variable mod protein filter and if
+   // current protein is not on the list, do not apply that particular variable mod.
+   // Any variable mod on the list will have pbVarModProteinFilter[?] = false unless
+   // the current protein matches a protein on the list.
+   if (g_staticParams.variableModParameters.mmapProteinLModsList.size() > 0)
+   {
+      char szProteinAccession[256];
+      sscanf(dbe.strName.c_str(), "%255s", szProteinAccession);
+      szProteinAccession[255] = '\0';
+
+      auto it = g_staticParams.variableModParameters.mmapProteinLModsList.begin();
+      while (it != g_staticParams.variableModParameters.mmapProteinLModsList.end())
+      {
+         int iWhichMod = it->first;
+
+         pbVarModProteinFilter[iWhichMod - 1] = false;  // do not apply this mod to this protein unless it's on the mmapProteinModsList
+
+         while (it != g_staticParams.variableModParameters.mmapProteinLModsList.end() && it->first == iWhichMod)
+         {
+            if (strstr(szProteinAccession, it->second.c_str()))
+            {
+               pbVarModProteinFilter[iWhichMod - 1] = true;
+            }
+            it++;
+         }
+      }
+   }
+
    // Quick clip n-term & PEFF variant check. Start summing amino acid mass at
    // the start variant position and work backwards.  If the mass is larger than
    // the max mass or the peptide length is longer than the max length before we
@@ -2124,6 +2161,17 @@ bool CometSearch::SearchForPeptides(struct sDBEntry dbe,
             // Otherwise, at this point, peptide mass is too big which means should be ok for varmod search.
             if (!g_staticParams.options.bCreateIndex && HasVariableMod(piVarModCounts, iStartPos, iEndPos, &dbe))
             {
+               // if variable mod protein filter applied, set residue mod count to 0 for the
+               // particular variable mod if current protein not on the protein filter list
+               if (g_staticParams.variableModParameters.mmapProteinLModsList.size() > 0)
+               {
+                  for (int i = 0; i < VMODS; ++i)
+                  {
+                     if (pbVarModProteinFilter[i] == false)
+                        piVarModCounts[i] = 0;
+                  }
+               }
+
                // VariableModSearch also includes looking at PEFF mods
                VariableModSearch(szProteinSeq, piVarModCounts, iStartPos, iEndPos, pbDuplFragment, &dbe);
             }

CometSearch/CometSearch.h

@@ -206,7 +206,7 @@ class CometSearch
                      int *piVarModSites,
                      struct sDBEntry *dbe);
    void VariableModSearch(char *szProteinSeq,
-                          int varModCounts[],
+                          int piVarModCounts[],
                           int iStartPos,
                           int iEndPos,
                           bool *pbDuplFragment,

CometSearch/CometSearchManager.cpp

@@ -624,7 +624,10 @@ bool CometSearchManager::InitializeStaticParams()
    {
       if (strData.length() > 0)
          strcpy(g_staticParams.szTxtFileExt, strData.c_str());
-   } 
+   }
+
+   if (GetParamValue("protein_modslist_file", strData))
+      g_staticParams.variableModParameters.sProteinLModsListFile = strData;
 
    if (GetParamValue("peff_obo", strData))
       strcpy(g_staticParams.peffInfo.szPeffOBO, strData.c_str());
@@ -805,8 +808,6 @@ bool CometSearchManager::InitializeStaticParams()
 
    GetParamValue("skip_researching", g_staticParams.options.bSkipAlreadyDone);
 
-// GetParamValue("skip_updatecheck", g_staticParams.options.bSkipUpdateCheck);
-
    GetParamValue("mango_search", g_staticParams.options.bMango);
 
    GetParamValue("scale_fragmentNL", g_staticParams.options.bScaleFragmentNL);
@@ -1982,13 +1983,59 @@ bool CometSearchManager::DoSearch()
    if (!g_staticParams.options.bOutputSqtStream) // && !g_staticParams.bIndexDb)
    {
       strOut = "\n Comet version \"" + g_sCometVersion + "\"\n\n";
-//      if (!g_staticParams.options.bSkipUpdateCheck)
-//       CometCheckForUpdates::CheckForUpdates(strOut.c_str());
 
       logout(strOut.c_str());
       fflush(stdout);
    }
 
+   FILE* fp;
+
+   // see if comet_varmod_proteins.txt is present, if so, read in the list of masses.
+   if ((fp = fopen(g_staticParams.variableModParameters.sProteinLModsListFile.c_str(), "r")) != NULL)
+   {
+      char szBuf[512];
+      vector<pair<int, string>> vpTmp;
+
+      printf(" Protein variable modifications filter:\n");
+
+      while (fgets(szBuf, 512, fp))
+      {
+         if (strlen(szBuf) > 3)
+         {
+            char szProtein[512];
+            int iWhichMod;
+
+            if (sscanf(szBuf, "%d %s", &iWhichMod, szProtein) == 2)
+            {
+               if (iWhichMod > 0 && iWhichMod <= VMODS)
+               {
+                  g_staticParams.variableModParameters.mmapProteinLModsList.insert({ iWhichMod, szProtein });
+               }
+            }
+         }
+      }
+      fclose(fp);
+
+      auto it = g_staticParams.variableModParameters.mmapProteinLModsList.begin();
+      while (it != g_staticParams.variableModParameters.mmapProteinLModsList.end())
+      {
+         int iWhichMod = it->first;
+         bool bFirst = true;
+
+         printf(" - variable_mod%02d: ", iWhichMod);
+         while (it != g_staticParams.variableModParameters.mmapProteinLModsList.end() && it->first == iWhichMod)
+         {
+            if (!bFirst)
+               printf(", ");
+            printf("%s", it->second.c_str());
+            it++;
+            bFirst = false;
+         }
+         printf("\n");
+      }
+      printf("\n");
+   }
+
    g_staticParams.precalcMasses.iMinus17 = BIN(g_staticParams.massUtility.dH2O);
    g_staticParams.precalcMasses.iMinus18 = BIN(g_staticParams.massUtility.dNH3);
    g_massRange.dMinMass = g_staticParams.options.dPeptideMassLow;