0
$\begingroup$

I am trying to calculate the substitution probabilities for amino acids (i.e. the probability of one amino acid given another). In theory, you should be able to derive them from a log-likelihood matrix (also called log-odds), such as Blosum62.

Does anyone know how to calculate probabilities from log-likelihoods or have a good modern source for a substitution probability matrix?

Improve this question
$\endgroup$
4
  • $\begingroup$ If I understand what you want I think you just need to align the sequences and then make some probability calculation (which will likely make no sense biologically but will work for you). If the sequences have rearrangements as I understand you said somewhere then it's more complicated than just a multiple-alignment and you need a tool that finds local alignments. Am I following you? $\endgroup$
    – juanjo75es
    Jan 22, 2021 at 23:55
  • $\begingroup$ I mean, if I understand well that's not going to be possible with any probabilistic method. You need to align first somehow. I guess something similar to what QUAST does. No idea if there is a tool for that. I said... $\endgroup$
    – juanjo75es
    Jan 23, 2021 at 0:02
  • $\begingroup$ Yes, the sequences are aligned first before comparing the amino acid probability of each position between the two groups. I was just wondering how to convert a Blosum62 log-likelihood score of "2" (for example) to a probability, because I can then use the set of probabilities to calculate a combined probability at every position for a set of sequences. $\endgroup$
    – chubs
    Jan 23, 2021 at 21:22
  • $\begingroup$ Note that you are using imprecise language. You don't specify probability of what. Anyway, I don't think I can help with that part. $\endgroup$
    – juanjo75es
    Jan 24, 2021 at 13:15

1 Answer 1

Reset to default
1
$\begingroup$

UPDATE2

@chubs wants to compare different clades. This is possible within recent releases of ConSurf:

consurf table 1

It is additionally quite possible to use ConSurf's default outputs to get per-position conservation scores that can be used for post-processing in any fashion, as well as prefab scripts that can be used with PyMol to make custom visualization:

For all cases, ConSurf creates the following outputs:

  1. The sequence and MSA colored by ConSurf conservation scores.

  2. A text file that summarizes for each position the normalized score calculated, the assigned grade (1-through-9), the reliability estimation (for the Bayesian method) and the amino acids/nucleotides observed in the respective MSA column.

  3. The sequences selected for the MSA and the MSA constructed (unless those files were uploaded by the user).

  4. A file with the frequency of each amino acid/nucleotide observed in each column of the MSA.

  5. The evolutionary tree, which was calculated by the server or uploaded by the user, together with the MSA are shown using the WASABI platform (Veidenberg et al., 2015). Moreover, using WASABI, a user can select a subtree containing a fraction of the homologous sequences and conduct a follow-up ConSurf analysis with these selected sequences. To refine a ConSurf analysis to a selected subtree, a user can choose any internal node on the phylogenetic tree and open a WASABI menu by using a right mouse click. Selecting the option 'run ConSurf on subtree' will issue a new window with the new ConSurf run for the selected subtree sequences (see example in Figure 3).

UPDATE

It seems that what @chubs really wants to do is map protein conservation onto a structure. There are various tools for doing this (here, here). ConSurf is the older and more established tool, using evolutionary conservation to rule out trivial apparent associations due to phylogenetic non-independence.

here is the somewhat related 3DPatch workflow, which uses profile HMM information as a measure of conservation:

3Dpatch workflow

Here is the output that they get:

3dpatch figure

Likewise ConSurf gets similar looking outputs:

consurf figure

ORIGINAL ANSWER

I think that more information is needed here to answer this question properly (what is your data? what is your ultimate goal? what are the things you've already tried?), but we can say a few things right away.

I would strongly suggest looking into substitution models, which are literally substitution probability matrices. Note that these are probabilities per unit time, which is to say that they are calibrated to e.g. evolutionary rates along a specific phylogenetic tree. I would therefore strongly suggest incorporating a phylogenetic model of some kind.

If you have e.g. a multiple sequence alignment, you can estimate a phylogenetic tree by one of the many tree estimation methods, which will give you a substitution probability matrix "for free" as this estimation is taken care of in the process of tree estimation. These will again be substitution probabilities per unit time; you will honestly not be able to get rid of the time-dependence here, as that's how evolution works.

I am a little behind the times, but some tools you might look into for tree estimation are PhyML and FastTree. If you only have sequences and haven't aligned them yet, you can use a tool like MUSCLE, ClustalO, or MAFFT.

BLOSUM and related models are derived from specific slices of protein sequence identity (evolutionary history) calibrated to certain scopes, it's unlikely that they are appropriate to your case, unless your data fits inside those rather narrow slices.

Improve this answer
$\endgroup$
15
  • $\begingroup$ I have a set of aligned sequences for which I want to generate a probability distribution of amino acids at each position. For example, at position 1, if all amino acids in the alignment are alanines, that position would hold a probability distribution for just alanine, which would be high for itself (alanine) and lower for other amino acids (just like the blosum62 row for alanine). If a position has different kinds of amino acids, its probability distribution would a mix, accordingly. $\endgroup$
    – chubs
    Jan 22, 2021 at 21:13
  • $\begingroup$ My ultimate goal is to have two seperate sets of sequences that can each be distilled to a probability distribution at every sequence position. I can then compare the two distributions (one from each set) to determine how similar that position is between the two groups. Sorry if this is confusing, but it's the best way I can think of for comparing two sets of alignments (e.g. mammalian sequences versus their arthropod counterparts). $\endgroup$
    – chubs
    Jan 22, 2021 at 21:18
  • $\begingroup$ @chubs is there a reason that you do not just use the sample distribution for each column of the alignment (for each taxonomic slice of alignment rows, if needed)? This is the normal methodology for e.g. generating logo plots, which I would recommend looking into: en.wikipedia.org/wiki/Sequence_logo If you are feeling very Bayesian you could pseudocount or try to incorporate a prior distribution based on the background frequency of amino acids in all positions of all the sequences, though that is getting a little over-elaborate IMO. $\endgroup$
    – Maximilian Press
    Jan 22, 2021 at 21:29
  • $\begingroup$ @chubs I'll note that methods such as PAML have specifically implemented phylogenetic modeling approaches for detecting amino acid patterns of the type that you seem to be interested in between clades, with the branch model allowing $\omega$ variation among branches: abacus.gene.ucl.ac.uk/software/paml.html $\endgroup$
    – Maximilian Press
    Jan 22, 2021 at 21:35
  • $\begingroup$ The reason I'm going about it this way is because I want to map the similarity of a given amino acid between two groups onto a three-dimensional atomic model of the protein. I've tried using sequence identity, but the results aren't quite what I want. In other words, one group with all glutamates at position 1 versus another group with all aspartates at the same position will give a 0% similarity when in fact this conservative substitution should give a higher score. I thought that by using probability distributions, I could calculate a more accurate similarity score to map onto the structure. $\endgroup$
    – chubs
    Jan 22, 2021 at 22:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.