Skip to content
Closed
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
5 changes: 3 additions & 2 deletions src/basic/basic.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.
#include "statistics.h"
#include "sequence.h"
#include "stats/standard_matrix.h"
#include "stats/matrices/accessors.h"
#include "util/log_stream.h"

const char* Const::version_string = "2.1.13";
Expand Down Expand Up @@ -278,7 +279,7 @@ Reduction::Reduction(const char* definition_string)
map_[letter] = i;
map8_[letter] = i;
map8b_[letter] = i;
freq_[i] += Stats::blosum62.background_freqs[letter];
freq_[i] += Stats::getBlosum62().background_freqs[letter];
}
for (double& f : freq_)
f = std::log(f);
Expand All @@ -298,4 +299,4 @@ std::string Reduction::decode_seed(const uint64_t seed, const size_t len) const
c /= size_;
}
return s;
}
}
32 changes: 32 additions & 0 deletions src/stats/matrices/accessors.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,32 @@
/****
DIAMOND protein aligner
Copyright (C) 2020 Max Planck Society for the Advancement of Science e.V.

Code developed by Benjamin Buchfink <benjamin.buchfink@tue.mpg.de>

This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
****/

#pragma once

namespace Stats {
const StandardMatrix& getBlosum45();
const StandardMatrix& getBlosum50();
const StandardMatrix& getBlosum62();
const StandardMatrix& getBlosum80();
const StandardMatrix& getBlosum90();
const StandardMatrix& getPam30();
const StandardMatrix& getPam70();
const StandardMatrix& getPam250();
}
10 changes: 7 additions & 3 deletions src/stats/matrices/blosum45.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix blosum45{
const StandardMatrix& getBlosum45() {

static const StandardMatrix blosum45_instance {

14, 2,

Expand Down Expand Up @@ -445,6 +447,8 @@ const StandardMatrix blosum45{
6.28722660e-01, 7.50000000e-01, 9.41772708e-01, 5.67250003e-01,
7.50000000e-01, 2.90000000e-01, 7.50000000e-01, 2.36320536e+00} }

};
};
return blosum45_instance;
}

}
}
10 changes: 7 additions & 3 deletions src/stats/matrices/blosum50.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix blosum50{
const StandardMatrix& getBlosum50() {

static const StandardMatrix blosum50_instance {

13, 2,

Expand Down Expand Up @@ -446,6 +448,8 @@ const StandardMatrix blosum50{
5.69808181e-01, 7.50000000e-01, 8.45421029e-01, 4.85898712e-01,
7.50000000e-01, 2.90000000e-01, 7.50000000e-01, 2.55759716e+00} }

};
};
return blosum50_instance;
}

}
}
11 changes: 8 additions & 3 deletions src/stats/matrices/blosum62.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix blosum62{
const StandardMatrix& getBlosum62() {

static const StandardMatrix blosum62_instance {

11, 1,

Expand Down Expand Up @@ -442,6 +444,9 @@ const StandardMatrix blosum62{
5.03239261e-01, 7.50000000e-01, 6.66952325e-01, 3.87839626e-01,
7.50000000e-01, 2.50000000e-01, 7.50000000e-01, 2.81516607e+00} }

};

}
};
return blosum62_instance;
}

}
8 changes: 6 additions & 2 deletions src/stats/matrices/blosum80.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix blosum80{
const StandardMatrix& getBlosum80() {

static const StandardMatrix blosum80_instance {

10, 1,

Expand Down Expand Up @@ -440,5 +442,7 @@ const StandardMatrix blosum80{
7.50000000e-01, 1.40000000e-01, 7.50000000e-01, 3.18242650e+00} }

};
return blosum80_instance;

}
}
}
10 changes: 7 additions & 3 deletions src/stats/matrices/blosum90.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix blosum90{
const StandardMatrix& getBlosum90() {

static const StandardMatrix blosum90_instance {

10, 1,

Expand Down Expand Up @@ -437,6 +439,8 @@ const StandardMatrix blosum90{
3.37037347e-01, 7.50000000e-01, 4.96615724e-01, 2.74820979e-01,
7.50000000e-01, 1.20000000e-01, 7.50000000e-01, 3.47015056e+00} }

};
};
return blosum90_instance;
}

}
}
10 changes: 7 additions & 3 deletions src/stats/matrices/pam250.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix pam250{
const StandardMatrix& getPam250() {

static const StandardMatrix pam250_instance {

14, 2,

Expand Down Expand Up @@ -445,6 +447,8 @@ const StandardMatrix pam250{
5.43022417e-01, 7.50000000e-01, 8.10366134e-01, 5.73996058e-01,
7.50000000e-01, 1.70000000e-01, 7.50000000e-01, 2.91088108e+00} }

};
};
return pam250_instance;
}

}
}
11 changes: 7 additions & 4 deletions src/stats/matrices/pam30.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix pam30{
const StandardMatrix& getPam30() {

static const StandardMatrix pam30_instance {

9, 1,

Expand Down Expand Up @@ -436,6 +438,7 @@ const StandardMatrix pam30{
9.17500230e-02, 7.50000000e-01, 1.02140077e-01, 1.07889016e-01,
7.50000000e-01, 3.00000000e-03, 7.50000000e-01, 6.85288369e+00} }

};

}
};
return pam30_instance;
}
}
8 changes: 6 additions & 2 deletions src/stats/matrices/pam70.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,9 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.

namespace Stats {

const StandardMatrix pam70{
const StandardMatrix& getPam70() {

static const StandardMatrix pam70_instance {

10, 1,

Expand Down Expand Up @@ -437,5 +439,7 @@ const StandardMatrix pam70{
7.50000000e-01, 2.00000000e-02, 7.50000000e-01, 5.58422761e+00} }

};
return pam70_instance;

}
}
}
7 changes: 4 additions & 3 deletions src/stats/score_matrix.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,7 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.
#include "score_matrix.h"
#include "basic/config.h"
#include "stats/cbs.h"
#include "stats/matrices/accessors.h"

using std::string;
using std::vector;
Expand Down Expand Up @@ -174,7 +175,7 @@ ScoreMatrix::ScoreMatrix(const string& matrix_file, int gap_open, int gap_extend
m[i][j] = score_array_[i * AMINO_ACID_COUNT + j];
p[i] = m[i];
}
const double* bg = Stats::blosum62.background_freqs.data();
const double* bg = Stats::getBlosum62().background_freqs.data();
try {
evaluer.initGapped(N, p, bg, bg, gap_open_, gap_extend_, gap_open_, gap_extend_, false, 0.01, 0.05, 120.0, 1024.0, 1);
}
Expand Down Expand Up @@ -238,12 +239,12 @@ void ScoreMatrix::init_background_scores()
for (size_t i = 0; i < 20; ++i) {
background_scores_[i] = 0;
for (size_t j = 0; j < 20; ++j)
background_scores_[i] += Stats::blosum62.background_freqs[j] * (*this)(i, j);
background_scores_[i] += Stats::getBlosum62().background_freqs[j] * (*this)(i, j);
}
}

double ScoreMatrix::bitscore(double raw_score) const
{
const double s = std::round(raw_score / scale_); // maintain compatibility with BLAST
return (lambda() * s - ln_k()) / LN_2;
}
}
7 changes: 4 additions & 3 deletions src/stats/stats.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -21,15 +21,16 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.
#include <algorithm>
#include <iterator>
#include "standard_matrix.h"
#include "matrices/accessors.h"

using std::string;
using std::max;
using std::min;

namespace Stats {

const std::map<std::string, const StandardMatrix&> StandardMatrix::matrices = { { "blosum45", blosum45 }, { "blosum62", blosum62 }, { "blosum50", blosum50 }, { "blosum80", blosum80 }, { "blosum90", blosum90 },
{ "pam250", pam250 }, { "pam30", pam30 }, { "pam70", pam70 } };
const std::map<std::string, const StandardMatrix&> StandardMatrix::matrices = { { "blosum45", getBlosum45() }, { "blosum62", getBlosum62() }, { "blosum50", getBlosum50() }, { "blosum80", getBlosum80() }, { "blosum90", getBlosum90() },
{ "pam250", getPam250() }, { "pam30", getPam30() }, { "pam70", getPam70() } };

const StandardMatrix& StandardMatrix::get(const std::string& name) {
string n;
Expand Down Expand Up @@ -59,4 +60,4 @@ double approx_id(Score raw_score, Loc range1, Loc range2) {
return min(max((double)raw_score / m * 16.56 + 11.41, 0.0), 100.0);
}

}
}
3 changes: 2 additions & 1 deletion src/tools/benchmark.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ along with this program. If not, see <http://www.gnu.org/licenses/>.
#include "dp/swipe/config.h"
#include "util/simd/dispatch.h"
#include "dp/score_vector_int16.h"
#include "stats/matrices/accessors.h"

void benchmark_io();

Expand Down Expand Up @@ -447,7 +448,7 @@ void matrix_adjust(const Sequence& s1, const Sequence& s2) {
high_resolution_clock::time_point t1 = high_resolution_clock::now();
vector<double> mat_final(TRUE_AA * TRUE_AA);
int iteration_count;
const double* joint_probs = (const double*)(Stats::blosum62.joint_probs);
const double* joint_probs = (const double*)(Stats::getBlosum62().joint_probs);
auto row_probs = Stats::composition(s1), col_probs = Stats::composition(s2);
config.cbs_err_tolerance = 0.0001;

Expand Down
Loading