Skip to content

Commit

Permalink
Site index starts from 0; introduce option --index-from-one to start …
Browse files Browse the repository at this point in the history 
…from 1
  • Loading branch information
@trongnhanuit
trongnhanuit committed Aug 16, 2023
1 parent a502340 commit 9224763
Show file tree
Hide file tree
Showing 5 changed files with 23 additions and 11 deletions.
2 changes: 1 addition & 1 deletion alignment/substitution.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@ Substitution::Substitution(const std::string& sub_str, Alignment* const aln, con
outError("Failed to parse the predefined mutation: '" + sub_str + "'. The new state is invalid.");

// Parse the position
position = convert_int(sub_str.substr(num_sites_per_state, input_length - (num_sites_per_state + num_sites_per_state)).c_str()) - 1;
position = convert_int(sub_str.substr(num_sites_per_state, input_length - (num_sites_per_state + num_sites_per_state)).c_str()) - Params::getInstance().site_starting_index;
if (aln->seq_type == SEQ_CODON)
position = position * ONE_THIRD;

Expand Down
14 changes: 7 additions & 7 deletions main/alisim.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -364,7 +364,7 @@ std::vector<std::pair<std::string,std::string>> readMutations(const std::string&
}
// if not found -> all characters are spaces -> invalid
else
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><old_state><site><new_state>");
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><list_of_mutations>");

// remove spaces at the ending of the line
pos = line.find_last_not_of(" ");
Expand All @@ -377,11 +377,11 @@ std::vector<std::pair<std::string,std::string>> readMutations(const std::string&
}
// if not found -> all characters are spaces -> invalid
else
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><old_state><site><new_state>");
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><list_of_mutations>");

// validate the input format <node_name><spaces><old_state><site><new_state>
// validate the input format <node_name><spaces><list_of_mutations>
if (line.length() < 5)
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><old_state><site><new_state>");
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><list_of_mutations>");

// replace \t by a space
line = regex_replace(line, std::regex("\t"), " ");
Expand All @@ -397,7 +397,7 @@ std::vector<std::pair<std::string,std::string>> readMutations(const std::string&
transform(node_name.begin(), node_name.end(), node_name.begin(), ::toupper);
}
else
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><old_state><site><new_state>");
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><list_of_mutations>");

// get mutation list at that node
std::string mutation_str = "";
Expand All @@ -408,7 +408,7 @@ std::vector<std::pair<std::string,std::string>> readMutations(const std::string&
if (line.length() - pos)
mutation_str = line.substr(pos, line.length() - pos);
else
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><old_state><site><new_state>");
outError("Invalid format '" + line + "'. Expected format should be <node_name><spaces><list_of_mutations>");

// record the mutations at that node
node_mutations.push_back(std::pair<std::string,std::string>(node_name, mutation_str));
Expand Down Expand Up @@ -646,7 +646,7 @@ void getLockedSites(Node* const node, Node* const dad, std::vector<bool>* const
// vailidate position
if (pos >= seq_length)
{
outWarning("Ignore a predefined mutation " + aln->convertStateBackStr(mut_it->getOldState()) + convertIntToString((aln->seq_type == SEQ_CODON ? pos * 3 : pos) + 1) + aln->convertStateBackStr(mut_it->getNewState()) + ". Position exceeds the sequence length " + convertIntToString(aln->seq_type == SEQ_CODON ? seq_length * 3 : seq_length));
outWarning("Ignore a predefined mutation " + aln->convertStateBackStr(mut_it->getOldState()) + convertIntToString((aln->seq_type == SEQ_CODON ? pos * 3 : pos) + Params::getInstance().site_starting_index) + aln->convertStateBackStr(mut_it->getNewState()) + ". Position exceeds the sequence length " + convertIntToString(aln->seq_type == SEQ_CODON ? seq_length * 3 : seq_length));
}
// mark the site as locked
else
Expand Down
6 changes: 3 additions & 3 deletions simulator/alisimulator.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1591,12 +1591,12 @@ void AliSimulator::handlePreMutations(const NeighborVec::iterator& it, int& pred
{
// show infor
if (verbose_mode >= VB_MED)
std::cout << "Apply a predefined mutation " << tree->aln->convertStateBackStr(mut_it->getOldState()) << convertIntToString(mut_it->getPosition() * num_sites_per_state + 1) << tree->aln->convertStateBackStr(mut_it->getNewState()) << std::endl;
std::cout << "Apply a predefined mutation " << tree->aln->convertStateBackStr(mut_it->getOldState()) << convertIntToString(mut_it->getPosition() * num_sites_per_state + params->site_starting_index) << tree->aln->convertStateBackStr(mut_it->getNewState()) << std::endl;

// if the old state is different from the one observed at the specified position -> show a warning
if ((*node_seq_chunk)[relative_pos] != mut_it->getOldState())
{
outWarning("State " + tree->aln->convertStateBackStr(mut_it->getOldState()) + " is not observed at position/site " + convertIntToString(mut_it->getPosition() * num_sites_per_state + 1) + ", we observe state " + tree->aln->convertStateBackStr((*node_seq_chunk)[relative_pos]) + " instead.");
outWarning("State " + tree->aln->convertStateBackStr(mut_it->getOldState()) + " is not observed at position/site " + convertIntToString(mut_it->getPosition() * num_sites_per_state + params->site_starting_index) + ", we observe state " + tree->aln->convertStateBackStr((*node_seq_chunk)[relative_pos]) + " instead.");
}

// apply the substitution
Expand All @@ -1608,7 +1608,7 @@ void AliSimulator::handlePreMutations(const NeighborVec::iterator& it, int& pred
// write log - for debugging only
/*ofstream myfile;
myfile.open ("thread_" + convertIntToString(omp_get_thread_num()) + ".txt", std::ostream::app);
myfile << tree->aln->convertStateBackStr(mut_it->getOldState()) << mut_it->getPosition() + 1 << tree->aln->convertStateBackStr(mut_it->getNewState())<< std::endl;
myfile << tree->aln->convertStateBackStr(mut_it->getOldState()) << mut_it->getPosition() + params->site_starting_index << tree->aln->convertStateBackStr(mut_it->getNewState())<< std::endl;
myfile.close();*/
}
}
Expand Down
7 changes: 7 additions & 0 deletions utils/tools.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1487,6 +1487,7 @@ void parseArg(int argc, char *argv[], Params &params) {
params.alignment_id = 0;
params.include_pre_mutations = false;
params.mutation_file = "";
params.site_starting_index = 0;

// store original params
for (cnt = 1; cnt < argc; cnt++) {
Expand Down Expand Up @@ -5267,6 +5268,12 @@ void parseArg(int argc, char *argv[], Params &params) {
continue;
}

if (strcmp(argv[cnt], "--index-from-one") == 0) {
params.site_starting_index = 1;

continue;
}

if (strcmp(argv[cnt], "--distribution") == 0) {
cnt++;
if (cnt >= argc || argv[cnt][0] == '-')
Expand Down
5 changes: 5 additions & 0 deletions utils/tools.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2601,6 +2601,11 @@ class Params {
* Mutation file that specifies pre-defined mutations occurs at nodes
*/
std::string mutation_file;

/**
* site starting index (for predefined mutations in AliSim)
*/
int site_starting_index;
};

/**
Expand Down

0 comments on commit 9224763

Please sign in to comment.