Loading...
Searching...
No Matches
freecontact.h
Go to the documentation of this file.
/* FreeContact - program to predict protein residue contacts from a sufficiently large multiple alignment
* Copyright (C) 2012-2013 Laszlo Kajan <lkajan@rostlab.org> Rost Lab, Technical University of Munich, Germany
*
* 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/>.
*/
#ifndef LIBFREECONTACT_H
#define LIBFREECONTACT_H
#include <errno.h>
#include <map>
#include <math.h>
#include <stdexcept>
#include <stdint.h>
#include <string.h>
#include <vector>
#ifdef __SSE2__
#include <x86intrin.h>
#endif // __SSE2__
#define LIBFREEC_API __attribute__ ((visibility ("default")))
#define LIBFREEC_LOCAL __attribute__ ((visibility ("hidden")))
namespace freecontact {
// http://www.fao.org/docrep/004/Y2775E/y2775e0e.htm
// Use 20 for gap.
static const uint8_t aamap_gapidx = 20;
static const uint8_t aamap[] = {
// A B C D E F G H I J K L M N
21, 0, 3, 4, 3, 6, 13, 7, 8, 9, 21, 11, 10, 12, 2,
// O P Q R S T U V W X Y Z
21, 14, 5, 1, 15, 16, 21, 19, 17, 21, 18, 6
};
static const uint8_t mapaa[] = {
//0 1 2 3 4 5 6 7 8 9
'A', 'R', 'N', 'D', 'C', 'Q', 'E', 'G', 'H', 'I',
'L', 'K', 'M', 'F', 'P', 'S', 'T', 'W', 'Y', 'V',
'-', 'X'
};
class LIBFREEC_API alilen_error : public std::runtime_error {
public:
alilen_error(const std::string& __arg) : std::runtime_error(__arg){}
};
class LIBFREEC_API icme_timeout_error : public std::runtime_error {
public:
icme_timeout_error(const std::string& __arg) : std::runtime_error(__arg){}
};
#ifndef __SSE2__
typedef long long __m128i __attribute__ ((__vector_size__ (16), __may_alias__));
#endif
class LIBFREEC_API ali_t : public std::vector<__m128i>
{
private:
#pragma GCC visibility push(hidden)
#pragma GCC visibility pop
public:
uint32_t seqcnt;
uint16_t alilen; // aligned sequence length, max 65535, reasonable
uint16_t alilen16; // size in 16-bytes
uint16_t alilenpad; // 16 * alilen16
ali_t( const uint16_t __alilen = 0 ) : std::vector<__m128i>(), seqcnt(0), alilen(__alilen), alilen16(( alilen-1 )/16 + 1), alilenpad( 16 * alilen16 )
{
reserve(1024*alilen16);
}
virtual ~ali_t(){}
inline uint8_t& operator()(uint32_t __k, uint16_t __ai)
{
return reinterpret_cast<uint8_t*>(this->_M_impl._M_start)[ __k*alilenpad + __ai ];
}
inline const uint8_t& operator()(uint32_t __k, uint16_t __ai) const
{
return reinterpret_cast<uint8_t*>(this->_M_impl._M_start)[ __k*alilenpad + __ai ];
}
static std::vector<uint8_t>
read_a_seq( const std::string& __l )
{
std::vector<uint8_t> ret; ret.reserve(__l.length());
for( std::string::const_iterator l_b = __l.begin(), l_e = __l.end(); l_b != l_e; ++l_b )
{
const std::string::value_type c = *l_b;
if( isalpha(c) )
ret.push_back( aamap[c & 0x1F] );
else if( c == '-' )
ret.push_back( 20 );
else break;
}
return ret;
}
ali_t& push(const std::vector<uint8_t>& __al);// throw (alilen_error);
inline ali_t& push(const std::string& __l) //throw (alilen_error)
{
return push(read_a_seq(__l));
}
};
struct LIBFREEC_API contact_t {
uint16_t i, j;
float score;
contact_t( uint16_t __i = 0, uint16_t __j = 0, float __score = 0 ) : i(__i), j(__j), score(__score) {}
inline bool operator<( const contact_t& __b ) const { return score < __b.score; } // to facilitate sorting by score
inline bool operator>( const contact_t& __b ) const { return score > __b.score; } // to facilitate sorting by score
};
struct parset_t {
double clustpc, density, gapth; uint16_t mincontsep;
double pseudocnt, pscnt_weight; bool estimate_ivcov; double shrink_lambda;
bool cov20, apply_gapth; double rho;
parset_t( double __clustpc = 0, double __density = 0, double __gapth = 0, uint16_t __mincontsep = 0,
double __pseudocnt = 0, double __pscnt_weight = 0, bool __estimate_ivcov = 0, double __shrink_lambda = 0,
bool __cov20 = 0, bool __apply_gapth = 0, double __rho = 0
) :
clustpc(__clustpc), density(__density), gapth(__gapth), mincontsep(__mincontsep),
pseudocnt(__pseudocnt), pscnt_weight(__pscnt_weight), estimate_ivcov(__estimate_ivcov), shrink_lambda(__shrink_lambda),
cov20(__cov20), apply_gapth(__apply_gapth), rho(__rho)
{}
};
const parset_t // clust conts gapt minc pseud pscw estim shrin cov20 agapth rho
//ps_evfold(0.70, 0.0, 1.0, 1, 0.0, 0.5, false, 0.0, true, false, -1 ), // original Matlab that expects gap pre-filtering
ps_evfold( 0.70, 0.0, 0.9, 1, 0.0, 0.5, false, 0.0, true, true, -1 ),
ps_psicov( 0.62, 0.03, 0.9, 5, 1.0, 0.0, true, 0.1, false, false, -1 ), // as published
ps_psicov_sd( 0.62, 0.0, 0.9, 5, 1.0, 0.0, true, 0.1, false, false, 0.001 ); // PSICOV sensible default (sd) for most cases
class LIBFREEC_API predictor {
public:
typedef std::map<std::string, std::vector<contact_t> > cont_res_t;
typedef double freq_t;
typedef float pairfreq_t;
typedef std::vector<freq_t> freq_vec_t;
typedef std::map<std::string, double> time_res_t;
public:
bool dbg;
predictor(bool __dbg = false) : dbg(__dbg) {}
virtual ~predictor(){}
void get_seq_weights(
freq_vec_t &__aliw, double &__wtot,
const ali_t& __ali, double __clustpc,
bool __veczw = true, int __num_threads = 0
);// throw (alilen_error, icme_timeout_error, std::range_error, std::runtime_error);
cont_res_t run(const ali_t& __ali, const freq_vec_t &__aliw, const double __wtot,
double __density, double __gapth, uint16_t __mincontsep,
double __pseudocnt, double __pscnt_weight, bool __estimate_ivcov, double __shrink_lambda,
bool __cov20, bool __apply_gapth, double __rho,
int __num_threads = 0, time_t __icme_timeout = 1800, time_res_t *__timing = NULL
);// throw (alilen_error, icme_timeout_error, std::range_error, std::runtime_error);
cont_res_t run(const ali_t& __ali, double __clustpc,
double __density, double __gapth, uint16_t __mincontsep,
double __pseudocnt, double __pscnt_weight, bool __estimate_ivcov, double __shrink_lambda,
bool __cov20, bool __apply_gapth, double __rho,
bool __veczw = true, int __num_threads = 0, time_t __icme_timeout = 1800, time_res_t *__timing = NULL
);// throw (alilen_error, icme_timeout_error, std::range_error, std::runtime_error);
inline
cont_res_t run(const ali_t& __ali,
const parset_t& __parset,
bool __veczw = true, int __num_threads = 0, time_t __icme_timeout = 1800, time_res_t *__timing = NULL
)// throw (alilen_error, icme_timeout_error, std::range_error, std::runtime_error)
{
return run(__ali, __parset.clustpc,
__parset.density, __parset.gapth, __parset.mincontsep,
__parset.pseudocnt, __parset.pscnt_weight, __parset.estimate_ivcov, __parset.shrink_lambda,
__parset.cov20, __parset.apply_gapth, __parset.rho,
__veczw, __num_threads, __icme_timeout, __timing);
}
};
} // namespace freecontact
#endif // LIBFREECONTACT_H
// vim:et:ts=4:ai:
freecontact::ali_t
Protein sequence alignment.
Definition freecontact.h:70
freecontact::ali_t::seqcnt
uint32_t seqcnt
Definition freecontact.h:75
freecontact::ali_t::push
ali_t & push(const std::string &__l)
Push a sequence to the alignment.
Definition freecontact.h:116
freecontact::ali_t::alilen
uint16_t alilen
Definition freecontact.h:76
freecontact::ali_t::push
ali_t & push(const std::vector< uint8_t > &__al)
Push a sequence to the alignment.
Definition freecontact.cpp:1001
freecontact::ali_t::~ali_t
virtual ~ali_t()
Definition freecontact.h:84
freecontact::ali_t::operator()
const uint8_t & operator()(uint32_t __k, uint16_t __ai) const
Definition freecontact.h:90
freecontact::ali_t::read_a_seq
static std::vector< uint8_t > read_a_seq(const std::string &__l)
Convert alignment string to amino acid codes with aamap.
Definition freecontact.h:97
freecontact::ali_t::alilen16
uint16_t alilen16
Definition freecontact.h:77
freecontact::ali_t::alilenpad
uint16_t alilenpad
Definition freecontact.h:78
freecontact::ali_t::ali_t
ali_t(const uint16_t __alilen=0)
Definition freecontact.h:80
freecontact::ali_t::operator()
uint8_t & operator()(uint32_t __k, uint16_t __ai)
Definition freecontact.h:86
freecontact::alilen_error::alilen_error
alilen_error(const std::string &__arg)
Definition freecontact.h:55
freecontact::icme_timeout_error::icme_timeout_error
icme_timeout_error(const std::string &__arg)
Definition freecontact.h:61
freecontact::predictor::cont_res_t
std::map< std::string, std::vector< contact_t > > cont_res_t
Definition freecontact.h:162
freecontact::predictor::dbg
bool dbg
Definition freecontact.h:168
freecontact::predictor::time_res_t
std::map< std::string, double > time_res_t
Definition freecontact.h:166
freecontact::predictor::~predictor
virtual ~predictor()
Definition freecontact.h:170
freecontact::predictor::pairfreq_t
float pairfreq_t
Definition freecontact.h:164
freecontact::predictor::run
cont_res_t run(const ali_t &__ali, const parset_t &__parset, bool __veczw=true, int __num_threads=0, time_t __icme_timeout=1800, time_res_t *__timing=NULL)
Predict residue contacts.
Definition freecontact.h:231
freecontact::predictor::predictor
predictor(bool __dbg=false)
Definition freecontact.h:169
freecontact::predictor::freq_vec_t
std::vector< freq_t > freq_vec_t
Definition freecontact.h:165
freecontact::predictor::freq_t
double freq_t
Definition freecontact.h:163
freecontact::predictor::run
cont_res_t run(const ali_t &__ali, const freq_vec_t &__aliw, const double __wtot, double __density, double __gapth, uint16_t __mincontsep, double __pseudocnt, double __pscnt_weight, bool __estimate_ivcov, double __shrink_lambda, bool __cov20, bool __apply_gapth, double __rho, int __num_threads=0, time_t __icme_timeout=1800, time_res_t *__timing=NULL)
Predict residue contacts.
Definition freecontact.cpp:434
LIBFREEC_API
#define LIBFREEC_API
Definition freecontact.h:31
freecontact
Definition freecontact.cpp:55
freecontact::ps_psicov_sd
const parset_t ps_psicov_sd(0.62, 0.0, 0.9, 5, 1.0, 0.0, true, 0.1, false, false, 0.001)
freecontact::aamap_gapidx
static const uint8_t aamap_gapidx
Definition freecontact.h:38
freecontact::__attribute__
long long __m128i __attribute__((__vector_size__(16), __may_alias__))
Definition freecontact.h:65
freecontact::ps_evfold
const parset_t ps_evfold(0.70, 0.0, 0.9, 1, 0.0, 0.5, false, 0.0, true, true, -1)
freecontact::ps_psicov
const parset_t ps_psicov(0.62, 0.03, 0.9, 5, 1.0, 0.0, true, 0.1, false, false, -1)
freecontact::aamap
static const uint8_t aamap[]
Definition freecontact.h:39
freecontact::mapaa
static const uint8_t mapaa[]
Definition freecontact.h:45
freecontact::contact_t::contact_t
contact_t(uint16_t __i=0, uint16_t __j=0, float __score=0)
Definition freecontact.h:129
freecontact::contact_t::operator<
bool operator<(const contact_t &__b) const
Definition freecontact.h:131
freecontact::contact_t::j
uint16_t j
Definition freecontact.h:126
freecontact::contact_t::score
float score
Contact score.
Definition freecontact.h:128
freecontact::contact_t::i
uint16_t i
Residue number, 0-based(!).
Definition freecontact.h:126
freecontact::contact_t::operator>
bool operator>(const contact_t &__b) const
Definition freecontact.h:132
freecontact::parset_t
Parameter set structure for predictor::run().
Definition freecontact.h:136
freecontact::parset_t::mincontsep
uint16_t mincontsep
Definition freecontact.h:137
freecontact::parset_t::cov20
bool cov20
Definition freecontact.h:139
freecontact::parset_t::shrink_lambda
double shrink_lambda
Definition freecontact.h:138
freecontact::parset_t::density
double density
Definition freecontact.h:137
freecontact::parset_t::rho
double rho
Definition freecontact.h:139
freecontact::parset_t::clustpc
double clustpc
Definition freecontact.h:137
freecontact::parset_t::estimate_ivcov
bool estimate_ivcov
Definition freecontact.h:138
freecontact::parset_t::gapth
double gapth
Definition freecontact.h:137
freecontact::parset_t::apply_gapth
bool apply_gapth
Definition freecontact.h:139
freecontact::parset_t::parset_t
parset_t(double __clustpc=0, double __density=0, double __gapth=0, uint16_t __mincontsep=0, double __pseudocnt=0, double __pscnt_weight=0, bool __estimate_ivcov=0, double __shrink_lambda=0, bool __cov20=0, bool __apply_gapth=0, double __rho=0)
Definition freecontact.h:141
freecontact::parset_t::pseudocnt
double pseudocnt
Definition freecontact.h:138
freecontact::parset_t::pscnt_weight
double pscnt_weight
Definition freecontact.h:138