Hola,
hoy he participado en un curso de introducción a SeqAn, una biblioteca Open Source escrita en lenguaje C++ para el análisis eficiente de secuencias biológicas, una tarea de la que ya habÃamos hablado en este blog. La versión original se publicó en BMC Bioinformatics y ha evolucionado hasta la actual SeqAn3, que puede obtenerse en https://github.com/seqan/seqan3 y require un compilador g++ con versión >= 7 y el estándar C++17 (g++-7 -std=c++17).
Para programar tus propias aplicaciones deberás tener algo de conocimiento de C++ moderno, en concreto de rangos, funciones lambda y la STL, pero afortunadamente el material del curso es muy ameno y fácil de seguir. Lo puedes encontrar en:
https://seqan.github.io/learning-resources/doc/biocpp/presentation
Os dejo un ejemplo que permite filtrar lecturas de secuencia de un fichero FASTQ que superen un cierto umbral de calidad y require leer la documentación de la API:
#include <string_view> #include <seqan3/std/algorithm> #include <seqan3/io/sequence_file/all.hpp> #include <seqan3/alphabet/quality/all.hpp> //ask user a PHRED quality cutoff value auto ask_quality() -> seqan3::phred42 { std::cout<<" Please type read quality cutoff [0,41]\n"; uint32_t cutoff{}; std::cin>>cutoff; return seqan3::assign_rank_to( cutoff, seqan3::phred42{}); } int main(int argc, const char * argv[]) { //note there's no user interface help nor proper arg checking std::string_view fastq_path = argv[1]; std::string_view fasta_path = argv[2]; // note file formats are deduced from file extensions seqan3::sequence_file_input fastq_file_in{fastq_path}; seqan3::sequence_file_output fasta_file_out{fasta_path}; // filter reads by quality seqan3::quality_alphabet qual_cutoff = ask_quality(); // this is where you need ti understand ranges & lambda functions, // plus a little bit about the SeqAnAPI for (const auto & [sequence,id,quality] : fastq_file_in) { if (std::ranges::all_of(quality, [qual_cutoff] (auto const & quality) { return quality >= qual_cutoff; })) fasta_file_out.emplace_back(sequence, id); } return 0; }
Si te interesa, el próximo curso gratuito será el 14 de septiembre de 16:00 a 19:00 CEST (más información en https://ogy.de/ts51). Hasta pronto,
Bruno
https://link.springer.com/book/10.1007/978-3-030-90990-1
ResponderEliminar