1 #ifndef _MUSSA_SEQUENCE_H_
2 #define _MUSSA_SEQUENCE_H_
3 // This file is part of the Mussa source distribution.
4 // http://mussa.caltech.edu/
5 // Contact author: Tristan De Buysscher, tristan@caltech.edu
7 // This program and all associated source code files are Copyright (C) 2005
8 // the California Institute of Technology, Pasadena, CA, 91125 USA. It is
9 // under the GNU Public License; please see the included LICENSE.txt
10 // file for more information, or contact Tristan directly.
13 // ----------------------------------------
14 // ---------- sequence.hh -----------
15 // ----------------------------------------
17 #include <boost/filesystem/path.hpp>
18 #include <boost/filesystem/fstream.hpp>
20 #include <boost/serialization/base_object.hpp>
21 #include <boost/serialization/list.hpp>
22 #include <boost/serialization/nvp.hpp>
23 #include <boost/serialization/string.hpp>
24 #include <boost/serialization/shared_ptr.hpp>
25 #include <boost/serialization/utility.hpp>
26 #include <boost/serialization/export.hpp>
28 #include <boost/shared_ptr.hpp>
35 // Sequence data class
37 //! Attach annotation information to a sequence track
41 annot(int begin, int end, std::string type, std::string name);
49 friend bool operator==(const annot& left, const annot& right);
51 // boost::serialization support
52 friend class boost::serialization::access;
53 template<class Archive>
54 void serialize(Archive& ar, const unsigned int /*version*/) {
55 ar & BOOST_SERIALIZATION_NVP(begin);
56 ar & BOOST_SERIALIZATION_NVP(end);
57 ar & BOOST_SERIALIZATION_NVP(type);
58 ar & BOOST_SERIALIZATION_NVP(name);
61 BOOST_CLASS_EXPORT(annot);
64 /* The way that motifs are found currently doesn't really
65 * indicate that the match was a reverse compliment
67 struct motif : public annot
71 motif() : annot(), sequence("") {};
72 //! this constructor is for when we're adding motifs to our annotations
73 motif(int begin, std::string motif);
76 // boost::serialization support
78 friend class boost::serialization::access;
79 template<class Archive>
80 void serialize(Archive& ar, const unsigned int /*version*/) {
81 ar & BOOST_SERIALIZATION_BASE_OBJECT_NVP(annot);
82 ar & BOOST_SERIALIZATION_NVP(sequence);
85 BOOST_CLASS_EXPORT(motif);
87 //! sequence track for mussa.
91 typedef std::string::difference_type difference_type;
92 typedef std::string::iterator iterator;
93 typedef std::string::const_iterator const_iterator;
94 typedef std::string::reference reference;
95 typedef std::string::const_reference const_reference;
96 typedef std::string::size_type size_type;
97 static const size_type npos = std::string::npos;
99 // some standard dna alphabets
100 // Include nl (\012), and cr (\015) to make sequence parsing eol
101 // convention independent.
103 static const std::string dna_alphabet;
104 static const std::string rna_alphabet;
105 //! this is the general iupac alphabet for nucleotides
106 static const std::string nucleic_iupac_alphabet;
107 //! the protein alphabet
108 static const std::string protein_alphabet;
112 Sequence(const char* seq);
113 Sequence(const std::string& seq);
114 Sequence(const Sequence& seq);
115 //! assignment to constant sequences
116 Sequence &operator=(const Sequence&);
118 // dangerous since they create large copies
119 //operator std::string();
120 //operator std::string() const;
122 friend std::ostream& operator<<(std::ostream&, const Sequence&);
123 friend bool operator<(const Sequence&, const Sequence&);
124 friend bool operator==(const Sequence&, const Sequence&);
125 const_reference operator[](size_type) const;
127 //! set sequence to a (sub)string containing nothing but AGCTN
128 void set_filtered_sequence(const std::string& seq,
129 std::string::size_type start=0,
130 std::string::size_type count=0);
132 //! retrive element at specific position
133 const_reference at(size_type n) const;
134 //! clear the sequence and its annotations
137 const_iterator begin() const;
139 const_iterator end() const;
140 //! is our sequence empty?
142 //! how many base pairs are there in our sequence
143 size_type size() const;
144 //! alias for size (used by string)
145 size_type length() const;
147 //! return a subsequence, copying over any appropriate annotation
148 Sequence subseq(int start=0, int count = std::string::npos) const;
149 //! return a reverse compliment
150 std::string rev_comp() const;
152 //! set sequence (filtered)
153 void set_sequence(const std::string &);
155 std::string get_sequence() const;
157 void set_species(const std::string &);
159 std::string get_species() const;
160 //! set the fasta header
161 void set_fasta_header(std::string header);
162 //! get the fasta header
163 std::string get_fasta_header() const;
164 //! get name (will return the first non-empty, of fasta_header, species)
165 std::string get_name() const;
167 //! load sequence AGCT from fasta file
168 //! \throw mussa_load_error
169 //! \throw sequence_empty_error
170 //! \throw sequence_empty_file_error
171 void load_fasta(const boost::filesystem::path file_path, int seq_num=1,
172 int start_index=0, int end_index=0);
173 //! load sequence from stream
174 //! \throw mussa_load_error
175 //! \throw sequence_empty_error
176 //! \throw sequence_empty_file_error
177 void load_fasta(std::iostream& file, int seq_num=1,
178 int start_index=0, int end_index=0);
179 //! load sequence annotations
180 //! \throws mussa_load_error
181 void load_annot(const boost::filesystem::path file_path, int start_index, int end_index);
182 //! load sequence annotations
183 //! \throws mussa_load_error
184 void load_annot(std::fstream& data_stream, int start_index, int end_index);
185 bool parse_annot(std::string data, int start_index=0, int end_index=0);
186 //! add an annotation to our list of annotations
187 void add_annotation(const annot& a);
188 const std::list<annot>& annotations() const;
189 const std::list<motif>& motifs() const;
191 //! add a motif to our list of motifs
192 //! \throws motif_normalize_error if there's something wrong with a_motif
193 void add_motif(const Sequence& a_motif);
194 //! clear our list of found motifs
196 //! search a sequence for a_motif
197 //! \throws motif_normalize_error if there's something wrong with a_motif
198 std::vector<int> find_motif(const std::string& a_motif) const;
199 //! search a sequence for a_motif
200 //! \throws motif_normalize_error if there's something wrong with a_motif
201 std::vector<int> find_motif(const Sequence& a_motif) const;
202 //! convert IUPAC symbols to upperase
203 //! \throws motif_normalize_error if there is an invalid symbol
204 static std::string motif_normalize(const std::string& a_motif);
206 //! annotate the current sequence with other sequences
207 void find_sequences(std::list<Sequence>::iterator start,
208 std::list<Sequence>::iterator end);
210 void save(boost::filesystem::fstream &save_file);
211 void load_museq(boost::filesystem::path load_file_path, int seq_num);
214 boost::shared_ptr<const std::string> seq;
218 std::list<annot> annots;
219 //! a seperate list for motifs since we're currently not saving them
220 std::list<motif> motif_list;
222 void motif_scan(std::string a_motif, std::vector<int> * motif_match_starts) const;
223 std::string rc_motif(std::string a_motif) const;
224 //! look for a string sequence type and and it to an annotation list
225 void add_string_annotation(std::string a_seq, std::string name);
227 // boost::serialization support
228 friend class boost::serialization::access;
229 template<class Archive>
230 void serialize(Archive& ar, const unsigned int /*version*/) {
231 ar & BOOST_SERIALIZATION_NVP(seq);
232 ar & BOOST_SERIALIZATION_NVP(header);
233 ar & BOOST_SERIALIZATION_NVP(species);
234 ar & BOOST_SERIALIZATION_NVP(annots);
235 ar & BOOST_SERIALIZATION_NVP(motif_list);
238 //BOOST_CLASS_EXPORT(Sequence);