Skip to content

Allow for single-state alignments and remove misleading warnings #73

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
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
42 changes: 23 additions & 19 deletions alignment/alignment.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -123,33 +123,37 @@ double chi2prob (int deg, double chi2)
} /* chi2prob */


int Alignment::checkAbsentStates(string msg) {
double *state_freq = new double[num_states];
computeStateFreq(state_freq);
string absent_states, rare_states;
int count = 0;
// Skip check for PoMo.
void Alignment::checkAbsentStates(string msg) {
// skip checking for PoMo
if (seq_type == SEQ_POMO)
return 0;
for (int i = 0; i < num_states; i++)
if (state_freq[i] == 0.0) {
return;
string absent_states, rare_states;
int absent_cnt = 0;
double *state_freqs = new double[num_states];
computeStateFreq(state_freqs);
for (int x = 0; x < num_states; ++x) {
if (state_freqs[x] == 0.0) {
if (!absent_states.empty())
absent_states += ", ";
absent_states += convertStateBackStr(i);
count++;
} else if (state_freq[i] <= Params::getInstance().min_state_freq) {
absent_states += convertStateBackStr(x);
absent_cnt++;
} else if (state_freqs[x] <= Params::getInstance().min_state_freq) {
if (!rare_states.empty())
rare_states += ", ";
rare_states += convertStateBackStr(i);
rare_states += convertStateBackStr(x);
}
if (count >= num_states-1 && Params::getInstance().fixed_branch_length != BRLEN_FIX)
outError("Only one state is observed in " + msg);
}
delete [] state_freqs;
if (absent_cnt == num_states)
outError("Only gaps observed in " + msg);
if (absent_cnt == num_states - 1)
outWarning("Only one state observed in " + msg);
if (absent_cnt > 0)
outWarning(convertIntToString(absent_cnt) + " states (see below) not observed in " + msg);
if (!absent_states.empty())
cout << "NOTE: State(s) " << absent_states << " not present in " << msg << " and thus removed from Markov process to prevent numerical problems" << endl;
outWarning("State(s) " + absent_states + " not present in " + msg + " and may cause numerical problems");
if (!rare_states.empty())
cout << "WARNING: States(s) " << rare_states << " rarely appear in " << msg << " and may cause numerical problems" << endl;
delete[] state_freq;
return count;
outWarning("State(s) " + rare_states + " rarely appear in " + msg + " and may cause numerical problems");
}

void Alignment::checkSeqName() {
Expand Down
5 changes: 2 additions & 3 deletions alignment/alignment.h
Original file line number Diff line number Diff line change
Expand Up @@ -484,11 +484,10 @@ class Alignment : public vector<Pattern>, public CharSet, public StateSpace {
int getMaxSeqNameLength();

/*
check if some state is absent, which may cause numerical issues
check if some states are absent, which may cause numerical issues
@param msg additional message into the warning
@return number of absent states in the alignment
*/
virtual int checkAbsentStates(string msg);
virtual void checkAbsentStates(string msg);

/**
check proper and undupplicated sequence names
Expand Down
8 changes: 3 additions & 5 deletions alignment/superalignment.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1148,11 +1148,9 @@ Alignment *SuperAlignment::removeIdenticalSeq(string not_remove, bool keep_two,
return aln;
}

int SuperAlignment::checkAbsentStates(string msg) {
int count = 0;
for (auto it = partitions.begin(); it != partitions.end(); ++it)
count += (*it)->checkAbsentStates("partition " + convertIntToString((it-partitions.begin())+1));
return count;
void SuperAlignment::checkAbsentStates(string msg) {
for (vector<Alignment*>::iterator it = partitions.begin(); it != partitions.end(); ++it)
(*it)->checkAbsentStates("partition " + convertIntToString((it-partitions.begin())+1));
}

/*
Expand Down
5 changes: 2 additions & 3 deletions alignment/superalignment.h
Original file line number Diff line number Diff line change
Expand Up @@ -171,11 +171,10 @@ class SuperAlignment : public Alignment


/*
check if some state is absent, which may cause numerical issues
check if some states are absent, which may cause numerical issues
@param msg additional message into the warning
@return number of absent states in the alignment
*/
virtual int checkAbsentStates(string msg);
virtual void checkAbsentStates(string msg);

/**
Quit if some sequences contain only gaps or missing data
Expand Down
18 changes: 2 additions & 16 deletions main/phyloanalysis.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -3403,22 +3403,8 @@ bool isTreeMixture(Params& params) {

void runTreeReconstruction(Params &params, IQTree* &iqtree) {

// string dist_file;
// params.startCPUTime = getCPUTime();
// params.start_real_time = getRealTime();

int absent_states = 0;
if (iqtree->isSuperTree()) {
PhyloSuperTree *stree = (PhyloSuperTree*)iqtree;
for (auto i = stree->begin(); i != stree->end(); i++)
absent_states += (*i)->aln->checkAbsentStates("partition " + (*i)->aln->name);
} else {
absent_states = iqtree->aln->checkAbsentStates("alignment");
}
if (absent_states > 0) {
cout << "NOTE: " << absent_states << " states (see above) are not present and thus removed from Markov process to prevent numerical problems" << endl;
}

iqtree->aln->checkAbsentStates("alignment");

// Make sure that no partial likelihood of IQ-TREE is initialized when PLL is used to save memory
if (params.pll) {
iqtree->deleteAllPartialLh();
Expand Down
Loading