OpenStructure
|
Namespaces | |
aaindex | |
mat | |
renumber | |
Data Structures | |
struct | AlignmentOpts |
struct | ContactPredictionScoreResult |
struct | ContactWeightMatrix |
struct | Distances |
class | DistanceMap |
class | ContextProfile |
class | ContextProfileDB |
class | InsDel |
struct | PairSubstWeightMatrix |
struct | RefMode |
class | SubstWeightMatrix |
class | VarianceMap |
class | Dist2Mean |
class | MeanlDDT |
Typedefs | |
typedef boost::shared_ptr< DistanceMap > | DistanceMapPtr |
typedef boost::shared_ptr< ContextProfileDB > | ContextProfileDBPtr |
typedef boost::shared_ptr< SubstWeightMatrix > | SubstWeightMatrixPtr |
typedef boost::shared_ptr< VarianceMap > | VarianceMapPtr |
typedef boost::shared_ptr< Dist2Mean > | Dist2MeanPtr |
typedef boost::shared_ptr< MeanlDDT > | MeanlDDTPtr |
Variables | |
char | RAW_CONTACT_WEIGHT_MATRIX_RES_LIST [21] ={'D', 'E', 'R', 'K', 'H', 'S', 'T', 'N', 'Q', 'G', 'P', 'Y', 'W', 'V', 'I', 'L', 'M', 'F', 'A', 'C','-'} |
Real | RAW_CONTACT_WEIGHT_MATRIX [21][21] |
char | RAW_PAIR_SUBST_WEIGHT_MATRIX_RES_LIST [20] ={'D', 'E', 'R', 'K', 'H', 'S', 'T', 'N', 'Q', 'G', 'P', 'Y', 'W', 'V', 'I', 'L', 'M', 'F', 'A', 'C'} |
Real | RAW_PAIR_SUBST_WEIGHT_MATRIX [20][20][20][20] |
typedef boost::shared_ptr<ContextProfileDB> ContextProfileDBPtr |
Definition at line 29 of file hmm_pseudo_counts.hh.
typedef boost::shared_ptr<Dist2Mean> Dist2MeanPtr |
Definition at line 41 of file variance_map.hh.
typedef boost::shared_ptr<DistanceMap> DistanceMapPtr |
Definition at line 148 of file distance_map.hh.
typedef boost::shared_ptr<MeanlDDT> MeanlDDTPtr |
Definition at line 42 of file variance_map.hh.
typedef boost::shared_ptr<SubstWeightMatrix> SubstWeightMatrixPtr |
Definition at line 36 of file subst_weight_matrix.hh.
typedef boost::shared_ptr<VarianceMap> VarianceMapPtr |
Definition at line 40 of file variance_map.hh.
void ost::seq::alg::AddAAPseudoCounts | ( | ost::seq::ProfileHandle & | profile, |
const ContextProfileDB & | db, | ||
Real | a = 0.9 , |
||
Real | b = 4.0 , |
||
Real | c = 1.0 |
||
) |
void ost::seq::alg::AddAAPseudoCounts | ( | ost::seq::ProfileHandle & | profile, |
Real | a = 1.0 , |
||
Real | b = 1.5 , |
||
Real | c = 1.0 |
||
) |
void ost::seq::alg::AddNullPseudoCounts | ( | ost::seq::ProfileHandle & | profile | ) |
void ost::seq::alg::AddTransitionPseudoCounts | ( | ost::seq::ProfileHandle & | profile, |
Real | gapb = 1.0 , |
||
Real | gabd = 0.15 , |
||
Real | gape = 1.0 |
||
) |
def ost.seq.alg.AlignmentFromChainView | ( | chain, | |
handle_seq_name = 'handle' , |
|||
view_seq_name = 'view' |
|||
) |
Creates and returns the sequence alignment of the given chain view to the chain handle. The alignment contains two sequences, the first containing all non-ligand peptide-linking residues, the second containing all non-ligand peptide-linking residues that are part of the view. :param chain: A valid chain :type chain: :class:`~ost.mol.ChainView` :param handle_seq_name: Name of the handle sequence in the output alignment :param view_seq_name: Name of the view sequence in the output alignment :returns: The alignment :rtype: :class:`~ost.seq.AlignmentHandle`
Definition at line 153 of file __init__.py.
def ost.seq.alg.AlignToSEQRES | ( | chain, | |
seqres, | |||
try_resnum_first = False , |
|||
validate = True |
|||
) |
Aligns the residues of chain to the SEQRES sequence, inserting gaps where needed. The function uses the connectivity of the protein backbone to find consecutive peptide fragments. These fragments are then aligned to the SEQRES sequence. All the non-ligand, peptide-linking residues of the chain must be listed in SEQRES. If there are any additional residues in the chain, the function raises a ValueError. :param chain: Source of the sequence :type chain: :class:`~ost.mol.ChainHandle` :param seqres: SEQRES sequence :type seqres: :class:`str` :param try_resnum_first: If set to True, this first builds an alignment using residue numbers and checks if the one-letter-codes match. If they all match, this alignment is used (and possibly validated). Otherwise, it displays a warning and falls back to the connectivity-based alignment. :type try_resnum_first: :class:`bool` :param validate: If set to True, the alignment is additionally checked by :func:`~ost.seq.alg.ValidateSEQRESAlignment` and raises a ValueError if the validation failed. :type validate: :class:`bool` :returns: The alignment of the residues in the chain and the SEQRES entries. :rtype: :class:`~ost.seq.AlignmentHandle`
Definition at line 60 of file __init__.py.
def ost.seq.alg.CalculateContactProbability | ( | cpred_res, | |
method | |||
) |
Calculate the probability of a predicted contact to be correct. This simply transforms the score associated with a prediction into a probability. :param cpred_res: A contact prediction result :param method: The method which was used for contact prediction. Should be one of "MI", "MIp", "MIpz", "cevoMI", "cevo" :type cpred_res: :class:`ost.seq.alg.ContactPredictionScoreResult` :type method: :class:`str`
Definition at line 223 of file __init__.py.
ContactPredictionScoreResult DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CalculateContactScore | ( | const AlignmentHandle & | aln, |
ContactWeightMatrix | w = LoadDefaultContactWeightMatrix() |
||
) |
ContactPredictionScoreResult DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CalculateContactSubstitutionScore | ( | const AlignmentHandle & | aln, |
int | ref_seq_index = 0 , |
||
PairSubstWeightMatrix | w = LoadDefaultPairSubstWeightMatrix() |
||
) |
ContactPredictionScoreResult DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CalculateMutualInformation | ( | const AlignmentHandle & | aln, |
ContactWeightMatrix | w = LoadConstantContactWeightMatrix() , |
||
bool | apc_correction = true , |
||
bool | zpx_transformation = true , |
||
float | small_number_correction = 0.05 |
||
) |
int DLLEXPORT_OST_SEQ_ALG ost::seq::alg::ClipAlignment | ( | AlignmentHandle & | aln, |
uint | n_seq_thresh = 2 , |
||
bool | set_offset = true , |
||
bool | remove_empty = true |
||
) |
Clips alignment so that first and last column have at least the desired number of structures.
aln | Multiple sequence alignment. Will be cut. |
n_seq_thresh | Minimal number of sequences desired. |
set_offset | Shall we update offsets for attached views? |
remove_empty | Shall we remove sequences with only gaps in cut aln? |
std::vector<Real> DLLEXPORT_OST_SEQ_ALG ost::seq::alg::Conservation | ( | const AlignmentHandle & | aln, |
bool | assign = true , |
||
const String & | prop_name = "cons" , |
||
bool | ignore_gap = false |
||
) |
Calculates conservation scores for each column in the alignment.
The conservation score is a value between 0 and 1. The bigger the number the more conserved the aligned residues are.
assign
If true, the conservation scores are assigned to attached residues. The name of the property can be changed with the prop_name parameter. Useful when coloring entities based on sequence conservation. prop_name
The property name for assigning the conservation to attached residues. Defaults to 'cons'. ignore_gap
If true, the dissimilarity between two gaps is increased to 6.0 instead of 0.5 as defined in the original version. Without this, a stretch where in the alignment there is only one sequence which is aligned to only gaps, is considered highly conserved (depending on the number of gap sequences).
Dist2MeanPtr DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CreateDist2Mean | ( | const DistanceMapPtr | dmap | ) |
dmap | Distance map as created with CreateDistanceMap. |
DistanceMapPtr DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CreateDistanceMap | ( | const seq::AlignmentHandle & | alignment | ) |
create distance map from a multiple sequence alignment.
The algorithm requires that the sequence alignment consists of at least two sequences. The sequence at index 0 serves as a frame of reference. All the other sequences must have an attached view and a properly set sequence offset.
For each of the attached views, the C-alpha distance pairs are extracted and mapped onto the corresponding C-alpha distances in the reference sequence.
MeanlDDTPtr DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CreateMeanlDDTHA | ( | const DistanceMapPtr | dmap | ) |
VarianceMapPtr DLLEXPORT_OST_SEQ_ALG ost::seq::alg::CreateVarianceMap | ( | const DistanceMapPtr | dmap, |
Real | sigma = 25 |
||
) |
dmap | Distance map as created with CreateDistanceMap. |
sigma | Used for weighting of variance measure (see Distances::GetWeightedStdDev) |
AlignmentList DLLIMPORT ost::seq::alg::GlobalAlign | ( | const ConstSequenceHandle & | s1, |
const ConstSequenceHandle & | s2, | ||
SubstWeightMatrixPtr & | subst, | ||
int | gap_open = -5 , |
||
int | gap_ext = -2 |
||
) |
Real ost::seq::alg::HMMScore | ( | const ost::seq::ProfileHandle & | profile_0, |
const ost::seq::ProfileHandle & | profile_1, | ||
const ost::seq::AlignmentHandle & | aln, | ||
int | seq_0_idx, | ||
int | seq_1_idx, | ||
Real | match_score_offset = -0.03 , |
||
Real | correl_score_weight = 0.1 , |
||
Real | del_start_penalty_factor = 0.6 , |
||
Real | del_extend_penalty_factor = 0.6 , |
||
Real | ins_start_penalty_factor = 0.6 , |
||
Real | ins_extend_penalty_factor = 0.6 |
||
) |
ContactWeightMatrix DLLEXPORT_OST_SEQ_ALG ost::seq::alg::LoadConstantContactWeightMatrix | ( | ) |
weight of 1 for all amino-acid pairs and 0 for gaps.
ContactWeightMatrix DLLEXPORT_OST_SEQ_ALG ost::seq::alg::LoadDefaultContactWeightMatrix | ( | ) |
statistical potential matrix containing interaction pseudo energies
PairSubstWeightMatrix DLLEXPORT_OST_SEQ_ALG ost::seq::alg::LoadDefaultPairSubstWeightMatrix | ( | ) |
AlignmentList DLLIMPORT ost::seq::alg::LocalAlign | ( | const ConstSequenceHandle & | s1, |
const ConstSequenceHandle & | s2, | ||
SubstWeightMatrixPtr & | subst, | ||
int | gap_open = -5 , |
||
int | gap_ext = -2 |
||
) |
AlignmentHandle DLLEXPORT_OST_SEQ_ALG ost::seq::alg::MergePairwiseAlignments | ( | const AlignmentList & | pairwise_alns, |
const ConstSequenceHandle & | ref_seq | ||
) |
merge a list of pairwise alignments into one multiple sequence alignment
All sequences in the pairwise sequence alignments are a realigned to the reference sequence. This is useful to merge the results of a BLAST or HMM database search into one multiple sequence alignment.
The method does not produce the optimal multiple sequence alignemnt for all the sequences.
pairwise_alns | is a list of AlignmentHandles, each containing two sequences |
ref_seq | is the reference sequence. The reference sequence must not contain any gaps. |
ost::seq::AlignmentList ost::seq::alg::ParaGlobalAlign | ( | const ost::seq::ConstSequenceHandle & | s1, |
const ost::seq::ConstSequenceHandle & | s2, | ||
ost::seq::alg::SubstWeightMatrixPtr & | subst, | ||
int | gap_open = -5 , |
||
int | gap_ext = -2 |
||
) |
ost::seq::AlignmentList ost::seq::alg::ParaLocalAlign | ( | const ost::seq::ConstSequenceHandle & | s1, |
const ost::seq::ConstSequenceHandle & | s2, | ||
ost::seq::alg::SubstWeightMatrixPtr & | subst, | ||
int | gap_open = -5 , |
||
int | gap_ext = -2 |
||
) |
bool ost::seq::alg::ParasailAvailable | ( | ) |
ost::seq::AlignmentList ost::seq::alg::ParaSemiGlobalAlign | ( | const ost::seq::ConstSequenceHandle & | s1, |
const ost::seq::ConstSequenceHandle & | s2, | ||
ost::seq::alg::SubstWeightMatrixPtr & | subst, | ||
int | gap_open = -5 , |
||
int | gap_ext = -2 |
||
) |
def ost.seq.alg.PredictContacts | ( | ali | ) |
Predicts contacts from a multiple sequence alignment using a combination of Mutual Information (*MI*) and the Contact Substitution Score (*CoEvoSc*). MI is calculated with the APC and small number corrections as well as with a transformation into Z-scores. The *CoEvoSc* is calculated using the default PairSubstWeightMatrix (see seq.alg.LoadDefaultPairSubstWeightMatrix). The final score for a pair of columns *(i,j)* of **ali** is obtained from: Sc(i,j)=MI(i,j)exp(CoEvoSc(i,j)) if *(i,j)* >=0 Sc(i,j)=MI(i,j)exp(1-CoEvoSc(i,j)) if *(i,j)* <0 :param ali: The multiple sequence alignment :type ali: :class:`~ost.seq.AlignmentHandle`
Definition at line 190 of file __init__.py.
AlignmentList DLLIMPORT ost::seq::alg::SemiGlobalAlign | ( | const ConstSequenceHandle & | s1, |
const ConstSequenceHandle & | s2, | ||
SubstWeightMatrixPtr & | subst, | ||
int | gap_open = -5 , |
||
int | gap_ext = -2 |
||
) |
Real DLLEXPORT_OST_SEQ_ALG ost::seq::alg::SequenceIdentity | ( | const AlignmentHandle & | aln, |
RefMode::Type | ref_mode = RefMode::LONGER_SEQUENCE , |
||
int | seq_a = 0 , |
||
int | seq_b = 1 |
||
) |
calculate sequence identity for two sequences in an alignment
ref_mode | influences the way the sequence identity is calculated. When set to LONGER_SEQUENCE, the sequence identity is calculated as the number of matches divided by the length of the longer sequence. If set to ALIGNMENT, the sequence identity is calculated as the number of matches divided by the number of aligned residues. |
seq_a | is the index of the first sequence |
seq_b | is the index of the second sequence |
aln | is the sequence alignment |
Real DLLEXPORT_OST_SEQ_ALG ost::seq::alg::SequenceSimilarity | ( | const AlignmentHandle & | aln, |
SubstWeightMatrixPtr | subst, | ||
bool | normalize = false , |
||
int | seq_a = 0 , |
||
int | seq_b = 1 |
||
) |
calculate sequence similarity for two sequences in an alignment
seq_a | is the index of the first sequence |
seq_b | is the index of the second sequence |
aln | is the sequence alignment |
std::vector<Real> DLLIMPORT ost::seq::alg::ShannonEntropy | ( | const AlignmentHandle & | aln, |
bool | ignore_gaps = true |
||
) |
calculates the Shannon entropy for each column in the alignment
def ost.seq.alg.ValidateSEQRESAlignment | ( | aln, | |
chain = None |
|||
) |
Checks if sequence in alignment has same connectivity as residues in chain. This looks for connected stretches in both the sequence and the chain and returns False if they don't match. This uses the connectivity of the protein backbone. :param aln: Alignment of two sequences with second one expected to map to residues in *chain*. :type aln: :class:`~ost.seq.AlignmentHandle` :param chain: Source of the sequence. :type chain: :class:`~ost.mol.ChainHandle` :returns: True if all residues (beside gapped ones) are connected, False otherwise.
Definition at line 4 of file __init__.py.
Real RAW_CONTACT_WEIGHT_MATRIX[21][21] |
Definition at line 5 of file default_contact_weight_matrix.hh.
char RAW_CONTACT_WEIGHT_MATRIX_RES_LIST[21] ={'D', 'E', 'R', 'K', 'H', 'S', 'T', 'N', 'Q', 'G', 'P', 'Y', 'W', 'V', 'I', 'L', 'M', 'F', 'A', 'C','-'} |
Definition at line 3 of file default_contact_weight_matrix.hh.
Real RAW_PAIR_SUBST_WEIGHT_MATRIX[20][20][20][20] |
Definition at line 5 of file default_pair_subst_weight_matrix.hh.
char RAW_PAIR_SUBST_WEIGHT_MATRIX_RES_LIST[20] ={'D', 'E', 'R', 'K', 'H', 'S', 'T', 'N', 'Q', 'G', 'P', 'Y', 'W', 'V', 'I', 'L', 'M', 'F', 'A', 'C'} |
Definition at line 3 of file default_pair_subst_weight_matrix.hh.