Merge pull request #262 from LadnerLab/prerelease

Make file and update for latest version

Author: jtladner

CMakeLists.txt

@@ -52,13 +52,24 @@ list( APPEND PepSIRF_LINK_LIBS
 )
 
 if(OpenMP_FOUND)
-  message( "OpenMP enabled" )
-  list( APPEND PepSIRF_LINK_LIBS OpenMP::OpenMP_CXX )
-  set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Xpreprocessor -fopenmp")
-  add_definitions( -DENABLE_OPENMP )
+  message("OpenMP enabled")
+  
+  if(APPLE)
+    # Get libomp filepath
+    execute_process(COMMAND brew --prefix libomp OUTPUT_VARIABLE BREW_PREFIX OUTPUT_STRIP_TRAILING_WHITESPACE)
+      
+    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Xpreprocessor -fopenmp -I${BREW_PREFIX}/include")
+    list(APPEND PepSIRF_LINK_LIBS "${BREW_PREFIX}/lib/libomp.dylib")
+  else()
+    list(APPEND PepSIRF_LINK_LIBS OpenMP::OpenMP_CXX)
+    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Xpreprocessor -fopenmp")
+  endif()
+
+  # Define OpenMP macro
+  add_definitions(-DENABLE_OPENMP)
 
 else()
-  message( "WARNING: OpenMP not found, parallelism disabled." )
+  message("WARNING: OpenMP not found, parallelism disabled.")
   set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wno-unknown-pragmas -Wno-unused-value")
 endif()
 

docs/5-changelog.md

@@ -10,10 +10,14 @@ permalink: /changelog/
 
 ## 1.7.0 | 2024-10-3
 
-<strong>Docker: added new feature (Issue #254).</strong> Added the ability to run PepSIRF as a Docker image and added a page for instructions.
+### Bug Fixes:
 
 <strong>CMakelists: bug fix (Issue #197).</strong> Resolved CMake not locating OpenMP on MacOS. Tutorial for fix added to installation page.
 
+## New Features:
+
+<strong>Docker: added new feature (Issue #254).</strong> Added the ability to run PepSIRF as a Docker image and added a page for instructions.
+
 <strong>Subjoin: added new feature (Issue #236).</strong> Added a functionality to the "-i" option in Subjoin to accept a regex pattern instead of a filename which contains sample/peptide names. The sample/peptide names used from the score matrix file will be filtered by whether they contain the regex pattern.
 
 <strong>Demux: added new feature (Issue #234).</strong> Added "--unmapped-reads-output" option to Demux, which writes all reads that have not been mapped to a sample/peptide to the specified filename.

extensions/e_k_bias.py

@@ -0,0 +1,69 @@
+#!/usr/bin/env python3
+
+import pandas as pd
+import argparse
+import fastatools as ft
+import numpy as np
+
+def main():
+
+	parser = argparse.ArgumentParser(formatter_class=argparse.ArgumentDefaultsHelpFormatter)
+ 
+	parser.add_argument('-i', '--fasta-file',  help='Directory with enriched petide files for input', required=True)
+	parser.add_argument('-o', '--output-file', default="e_k_bias_out.tsv", help='Name of .tsv to output file with AA bias data')
+	
+	args = parser.parse_args()
+
+	# get proportion of e's and k's for each peptide
+	e_k_props = get_e_k_props(args.fasta_file)
+
+	# get percentiles
+	percentile_dict = get_percentiles(e_k_props)
+
+	# create output df
+	out_data = [(name, round(prop, 3), round(percentile_dict[name], 2)) for name, prop in e_k_props.items()]
+
+	pd.DataFrame(out_data, columns=["CodeName", "e_k_Prop", "e_k_Percentile"]).to_csv(args.output_file, index=False, sep='\t')
+
+
+# get proportion of e's and k's for each peptide
+def get_e_k_props(fasta_file)->dict:
+	e_k_props = dict()
+
+	# get props for peptide file
+	fasta_dict = ft.read_fasta_dict(fasta_file)
+
+	# iterate through each sequence
+	for name, seq in fasta_dict.items():
+		e_k_count = 0
+
+		# loop through each AA, get count of e and k
+		for aa in seq:
+			if aa.lower() == 'e' or aa.lower() =='k':
+				e_k_count += 1
+
+		# add proportion to dict
+		e_k_props[name] = (e_k_count) / len(seq)
+
+	return e_k_props
+
+# get percentile of each peptide using its e and k proportion
+def get_percentiles(e_k_props)->dict:
+	# Calculate percentile of each peptide
+	names = list(e_k_props.keys())
+	all_props = np.array(list(e_k_props.values()))
+
+	# Get unique values and array for mapping each original value to its corresponding index in the unique array
+	unique_props, inverse_indices = np.unique(all_props, return_inverse=True)
+
+	# Calculate percentiles based on the unique props
+	percentile_ranks = np.linspace(0, 100, len(unique_props))
+
+	# Map sorted values back to the original corresponding name
+	percentile_dict = {names[i]: percentile_ranks[inverse_indices[i]] for i in range(len(names))}
+
+	return percentile_dict
+
+
+if __name__ == "__main__":
+	main()

extensions/findEpitopes.py

@@ -223,6 +223,7 @@ def read_check_align_file(directory):
 
     # Construct the full file path
     filepath = os.path.join(directory, 'checkAlignLength.out')
+
     #filepath = os.path.join('checkAlignLength.out')
 
     # Read the file content
@@ -253,6 +254,7 @@ def process_files_probes(probes_dict, directory_path):
         aligned_probes_file = filename.replace('.fasta', '_probesAligned.txt')
         aligned_probes_path = os.path.join(directory_path, aligned_probes_file)
 
+
         result[filename], file_2_pep_pos_dict[filename] = process_file_probes(data, aligned_probes_path)
 
     return result, file_2_pep_pos_dict

include/modules/core/matrix.h

@@ -85,7 +85,7 @@ class matrix
         bool end()
         {
             return current_idx >=
-                 this->matr.access_to_1d( this->matr.N, this->matr.M );
+                 this->matr->access_to_1d( this->matr->N, this->matr->M );
         }
 
         /**
@@ -126,7 +126,7 @@ class matrix
          **/
         mutable_iterator& next_row()
             {
-                this += this->matr.M;
+                this += this->matr->M;
                 return *this;
             }
 
@@ -139,7 +139,7 @@ class matrix
          **/
         mutable_iterator &next_col()
             {
-                this += this->matr.N;
+                this += this->matr->N;
                 return *this;
             }
 
@@ -240,7 +240,7 @@ class matrix
         bool end()
         {
             return current_idx >=
-                 this->matr.access_to_1d( this->matr.N, this->matr.M );
+                 this->matr->access_to_1d( this->matr->N, this->matr->M );
         }
 
         /**
@@ -281,7 +281,7 @@ class matrix
          **/
         iterator& next_row()
             {
-                this += this->matr.M;
+                this += this->matr->M;
                 return *this;
             }
 
@@ -294,7 +294,7 @@ class matrix
          **/
         iterator &next_col()
             {
-                this += this->matr.N;
+                this += this->matr->N;
                 return *this;
             }
 

include/modules/core/pepsirf_version.h

@@ -1,6 +1,6 @@
 #ifndef PEPSIRF_VERSION_HH_INCLUDED
 #define PEPSIRF_VERSION_HH_INCLUDED
 
-#define PEPSIRF_VERSION "1.6.0"
+#define PEPSIRF_VERSION "1.7.0"
 
 #endif

include/modules/core/sequence.h

@@ -72,6 +72,14 @@ class sequence
      **/
     bool operator==( const sequence& s ) const;
 
+    /**
+     * Less than operator override. 
+     * For two sequences a and b we say a < b iff 
+     * a.seq < b.seq. 
+     * @param s Sequence to compare against.
+     **/
+    bool operator<( const sequence& s ) const;
+
 
     /**
      * Get the length of the 'seq' member of this class.

include/modules/demux/module_demux.h

@@ -109,6 +109,22 @@ class module_demux : public module
         std::vector<sequence> &lib_seqs
     );
 
+    /**
+     * Outputs truncated sequence info for unqiue sequences, non-unqiue 
+     * sequences, and the new fasta-formatted file
+     * @param seq_length Sequence length specified with the "--seq" option
+     * @param lib_seqs Value of library sequences received from file
+     *        specified by the "--library"
+     * @param library_fname Library file name, specified by "--library"
+     * @param trunc_info_outdir Specified output directory with the 
+     *        "--trunc_info_outdir" option
+     **/
+    void output_trunc_info(
+        std::size_t seq_length,
+        std::vector<sequence> lib_seqs,
+        std::string library_fname,
+        std::string trunc_info_outdir
+    );
 
     /**
      * Method to zero a vector of size_t elements.
@@ -380,11 +396,32 @@ class module_demux : public module
      * Creates a single output fastq file containing all of the reads that have not been mapped to a sample/peptide
      * @param filename file to output to
      * @param samp_map fastaq output map
-     * @reads_dup vector of all reads
+     * @param reads_dup vector of all reads
      **/
     void create_unmapped_reads_file( std::string filename, 
                             std::map<std::string, std::vector<fastq_sequence>> samp_map, std::vector<fastq_sequence> reads_dup );
 
+    /**
+     * Outputs a vector to a file with a delemiter between entries
+     * @param vector vector to be outputted, make sure the type of its values can use the << operator for output
+     * @param filename file to output to
+     * @param delimiter delimiter to output between values
+     **/
+    template <typename T>
+    void simple_vector_out(std::vector<T> vector, std::string output_file, std::string delimiter)
+    {
+        std::ofstream outfile(output_file, std::ios::out);
+        for (unsigned int i=0; i<vector.size(); i++)
+        {
+            outfile << vector[i];
+            if( i < vector.size()-1 )
+            {
+                outfile << delimiter;
+            }
+        }
+        outfile.close();
+    }
+
 };
 
 #endif /* MODULE_DEMUX_HH_INCLUDED */

include/modules/demux/options_demux.h

@@ -36,9 +36,9 @@ class options_demux: public options
     int num_indexes;
     std::string replicate_info_fname;
     std::string unmapped_reads_fname;
-
     bool translation_aggregation;
     std::string fastq_out;
+    std::string trunc_info_outdir;
 
     /**
      * The number of fastq records to read per loop. A higher value here will result in higher memory usage by the program.

src/modules/core/sequence.cpp

@@ -53,3 +53,8 @@ bool sequence::operator==( const sequence& s ) const
 {
     return !s.seq.compare( seq );
 }
+
+bool sequence::operator<( const sequence& s ) const
+{
+    return seq < s.seq;
+}