diff --git a/modules/seq/base/doc/seq.rst b/modules/seq/base/doc/seq.rst
index cecd2fd1bc451497034bd8e3d415c52d1b87500f..da94a428bcd6f589c9cab4ae2c084c597de13e10 100644
--- a/modules/seq/base/doc/seq.rst
+++ b/modules/seq/base/doc/seq.rst
@@ -503,6 +503,10 @@ single columns containing amino acid frequencies and transition probabilities.
     :param col:  Column to add
     :type col:  :class:`HMMColumn`
 
+  .. attribute:: sequence
+
+    Sequence of the columns
+
   .. attribute:: columns
 
     Iterable columns of the HMM
diff --git a/modules/seq/base/pymod/export_hmm.cc b/modules/seq/base/pymod/export_hmm.cc
index 944892f80237974f21463e7939d47aea5c67688b..fea380dcab1f3fb56d3077722d125f5b9020a55f 100644
--- a/modules/seq/base/pymod/export_hmm.cc
+++ b/modules/seq/base/pymod/export_hmm.cc
@@ -44,5 +44,6 @@ void export_hmm()
                   make_function(&HMM::GetColumns,
                   return_value_policy<copy_const_reference>()))
     .add_property("avg_entropy", &HMM::GetAverageEntropy)
+    .add_property("sequence",&HMM::GetSequence)
   ;
 }
diff --git a/modules/seq/base/src/hmm.cc b/modules/seq/base/src/hmm.cc
index 07555a3c5bccac942a7e4912d6a92bab3b4a1dbb..e53b062731f668fd06bfcc5b17d41bad86c2fa93 100644
--- a/modules/seq/base/src/hmm.cc
+++ b/modules/seq/base/src/hmm.cc
@@ -149,6 +149,16 @@ HMMPtr HMM::Load(const std::string& filename) {
   return hmm;
 }
 
+String HMM::GetSequence() const{
+
+  std::stringstream ss;
+  for(HMMColumnList::const_iterator i = this->columns_begin(); 
+      i != this->columns_end(); ++i){
+    ss << i->GetOneLetterCode();
+  }
+  return ss.str();
+}
+
 
 Real HMM::GetAverageEntropy() const {
   Real n_eff=0.0;
diff --git a/modules/seq/base/src/hmm.hh b/modules/seq/base/src/hmm.hh
index 9f1871827a730a57b62e51d034f162b17f559b1a..0dcaf36866e316ca1bd53b13d44406252fc21845 100644
--- a/modules/seq/base/src/hmm.hh
+++ b/modules/seq/base/src/hmm.hh
@@ -17,26 +17,36 @@ typedef enum {
 
 class HMMColumn {
  public:
+
   HMMColumn() : n_eff_(0.0), n_eff_ins_(0.0), n_eff_del_(0.0) { 
     memset(freq_, 0, sizeof(Real)*20);
     memset(trans_, 0, sizeof(Real)*9);
   }
+
   HMMColumn(const HMMColumn& rhs): olc_(rhs.olc_), n_eff_(rhs.n_eff_), 
     n_eff_ins_(rhs.n_eff_ins_), n_eff_del_(rhs.n_eff_del_) {
     memcpy(freq_, rhs.freq_, sizeof(Real)*20);
     memcpy(trans_, rhs.trans_, sizeof(Real)*9);
   }
+
   Real GetTransitionFreq(HMMState from, HMMState to) const {
     return trans_[from][to];
   }
+
   void SetTransitionFreq(HMMState from, HMMState to, Real freq){
     trans_[from][to] = freq;
   }
+
   void SetNEff(Real val) { n_eff_ = val; }
+
   void SetNEffIns(Real val) { n_eff_ins_ = val; }
+
   void SetNEffDel(Real val) { n_eff_del_ = val; }
+
   Real GetNEff() const { return n_eff_; }
+
   Real GetNEffIns() const { return n_eff_ins_; }
+
   Real GetNEffDel() const { return n_eff_del_; }
 
   Real GetFreq(char ch) const {
@@ -169,11 +179,17 @@ typedef std::vector<HMMColumn> HMMColumnList;
 class HMM { 
  public:
   HMM() {}
+
   static HMMPtr Load(const std::string& filename);
+
   const std::vector<HMMColumn>& GetColumns() const { return columns_; }
+
   const HMMColumn& GetNullModel() const { return null_model_; } 
+
   void SetNullModel(const HMMColumn& null_model) { null_model_ = null_model; }
 
+  String GetSequence() const;
+
   //some functions to make it behave like a vector
 
   size_t size() const { return columns_.size(); }
@@ -200,6 +216,7 @@ class HMM {
   HMMColumnList::iterator columns_end() { return columns_.end(); }
   HMMColumnList::const_iterator columns_begin() const { return columns_.begin(); }
   HMMColumnList::iterator columns_begin() { return columns_.begin(); }
+
   Real GetAverageEntropy() const;
 
   friend std::ofstream& operator<<(std::ofstream& os, HMM& hmm){