Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
68 changes: 68 additions & 0 deletions src/seqTools.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,68 @@
//
// seqTools.cpp
// StochHMM
//
// Created by Paul Lott on 5/18/12.
// Copyright (c) 2012 University of California, Davis. All rights reserved.
//

#include "seqTools.h"

namespace StochHMM{


sequence shuffle(sequence *seq){
sequence shuffled_seq(ALPHA_NUM);
std::string strSeq=seq->getUndigitized();
for(size_t i=0;i<strSeq.size();++i){
size_t k = rand() % (strSeq.size());
if (k!=i){
iter_swap(strSeq.begin()+i, strSeq.begin()+k);
}
}

shuffled_seq.setSeq(strSeq, seq->getTrack());

return shuffled_seq;
}

sequence random_sequence(std::vector<double>& freq, size_t length, track* tr){
sequence random_seq(ALPHA_NUM);

if (tr==NULL){
std::cerr << "Track is not defined" << std::endl;
return random_seq;
}

size_t alphaSize=tr->getAlphaSize();
size_t freqSize=freq.size();
if (alphaSize!=freqSize){
std::cerr << "Frequency distribution size and Alphabet size must be the same." << std::endl;
return random_seq;
}

//Create CDF of frequency distribution
std::vector<std::pair<double,std::string> > cdf;
double sum = 0.0;
for(size_t i=0;i<freqSize;++i){
sum+=freq[i];
std::pair<double,std::string> val (sum, tr->getAlpha(i));
cdf.push_back(val);
}

//Generate random sequence
std::string random_string;
for(size_t j=0;j<length;++j){
double val = ((double)rand()/((double)(RAND_MAX)+(double)(1))); //Generate random
for (size_t m=0;m<freqSize;++m){ //Check to see which value is
if (cdf[m].first>=val){
random_string+=cdf[m].second;
break;
}
}
}

random_seq.setSeq(random_string, tr);
return random_seq;
}
}
39 changes: 39 additions & 0 deletions src/seqTools.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,39 @@
//
// seqTools.h
// StochHMM
//
// Created by Paul Lott on 5/18/12.
// Copyright (c) 2012 University of California, Davis. All rights reserved.
//

#ifndef StochHMM_seqTools_h
#define StochHMM_seqTools_h
#include <iostream>
#include <string>
#include <algorithm>
#include "sequence.h"
#include "hmm.h"
namespace StochHMM{
sequence shuffle(sequence* seq);

sequence random_sequence(std::vector<double> const& freq, size_t , track*);
sequence random_sequence(emm*);

sequence reverseComplement();
sequence translate();


void motifScoring();

void markov_length_distribution(model*);

void markov_generate_sequence(model*);

//Put in PWM scoring with options for set threshold or determine threshold through shuffling and calculation of FDR(give an FDR threshold)


}



#endif