Alineamiento de secuencias y matrices de sustitución

Muchos de vosotros ya sabréis de la importancia de los alineamientos de secuencias en bioinformática, puesto que tienen aplicaciones en cosas tan diferentes como la predicción de genes en genomas o la predicción de estructura secundaria en proteínas.

Los algoritmos de alineamiento de secuencias, como BLAST, usan matrices de sustitución para dirigir el proceso de alineamiento. Para proteínas, una de las matrices más conocidas es BLOSUM62.


Tabla 2.1: BLOSUM62, muestra estimaciones del valor evolutivo del cambio de un aminoácido por otro en una proteína. Por ejemplo, una mutación de alanina (A) a aspártico (D) está penalizada con un -2.
A R N D C Q E G H I L K M F P S T W Y V B Z X *
A 4 -1 -2 -2 0 -1 -1 0 -2 -1 -1 -1 -1 -2 -1 1 0 -3 -2 0 -2 -1 0 -4
R -1 5 0 -2 -3 1 0 -2 0 -3 -2 2 -1 -3 -2 -1 -1 -3 -2 -3 -1 0 -1 -4
N -2 0 6 1 -3 0 0 0 1 -3 -3 0 -2 -3 -2 1 0 -4 -2 -3 3 0 -1 -4
D -2 -2 1 6 -3 0 2 -1 -1 -3 -4 -1 -3 -3 -1 0 -1 -4 -3 -3 4 1 -1 -4
C 0 -3 -3 -3 9 -3 -4 -3 -3 -1 -1 -3 -1 -2 -3 -1 -1 -2 -2 -1 -3 -3 -2 -4
Q -1 1 0 0 -3 5 2 -2 0 -3 -2 1 0 -3 -1 0 -1 -2 -1 -2 0 3 -1 -4
E -1 0 0 2 -4 2 5 -2 0 -3 -3 1 -2 -3 -1 0 -1 -3 -2 -2 1 4 -1 -4
G 0 -2 0 -1 -3 -2 -2 6 -2 -4 -4 -2 -3 -3 -2 0 -2 -2 -3 -3 -1 -2 -1 -4
H -2 0 1 -1 -3 0 0 -2 8 -3 -3 -1 -2 -1 -2 -1 -2 -2 2 -3 0 0 -1 -4
I -1 -3 -3 -3 -1 -3 -3 -4 -3 4 2 -3 1 0 -3 -2 -1 -3 -1 3 -3 -3 -1 -4
L -1 -2 -3 -4 -1 -2 -3 -4 -3 2 4 -2 2 0 -3 -2 -1 -2 -1 1 -4 -3 -1 -4
K -1 2 0 -1 -3 1 1 -2 -1 -3 -2 5 -1 -3 -1 0 -1 -3 -2 -2 0 1 -1 -4
M -1 -1 -2 -3 -1 0 -2 -3 -2 1 2 -1 5 0 -2 -1 -1 -1 -1 1 -3 -1 -1 -4
F -2 -3 -3 -3 -2 -3 -3 -3 -1 0 0 -3 0 6 -4 -2 -2 1 3 -1 -3 -3 -1 -4
P -1 -2 -2 -1 -3 -1 -1 -2 -2 -3 -3 -1 -2 -4 7 -1 -1 -4 -3 -2 -2 -1 -2 -4
S 1 -1 1 0 -1 0 0 0 -1 -2 -2 0 -1 -2 -1 4 1 -3 -2 -2 0 0 0 -4
T 0 -1 0 -1 -1 -1 -1 -2 -2 -1 -1 -1 -1 -2 -1 1 5 -2 -2 0 -1 -1 0 -4
W -3 -3 -4 -4 -2 -2 -3 -2 -2 -3 -2 -3 -1 1 -4 -3 -2 11 2 -3 -4 -3 -2 -4
Y -2 -2 -2 -3 -2 -1 -2 -3 2 -1 -1 -2 -1 3 -3 -2 -2 2 7 -1 -3 -2 -1 -4
V 0 -3 -3 -3 -1 -2 -2 -3 -3 3 1 -2 1 -1 -2 -2 0 -3 -1 4 -3 -2 -1 -4
B -2 -1 3 4 -3 0 1 -1 0 -3 -4 0 -3 -3 -2 0 -1 -4 -3 -3 4 1 -1 -4
Z -1 0 0 1 -3 3 4 -2 0 -3 -3 1 -1 -3 -1 0 -1 -3 -2 -2 1 4 -1 -4
X 0 -1 -1 -1 -2 -1 -1 -1 -1 -1 -1 -1 -1 -1 -2 0 0 -2 -1 -1 -1 -1 -1 -4
* -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 1


En este ejemplo vamos a construir una matriz de sustitución de aminoácidos, como BLOSUM62, pero basada en el código genético. La hipótesis que vamos a manejar es que los aminoácidos cuyos codones se parecen más son más fácilmente intercambiables, olvidándonos por un momento de sus propiedades bioquímicas.

 
#!/usr/bin/perl -w 
# Ejercicio escrito por Bruno Contreras

# Programa que calcula una matriz de sustitucion de aminoacidos en base al parecido de 
# sus codones de tRNA

use strict; 

my %tripletes = (
'F' => ['TTT','TTC'],
'L' => ['TTA','TTG','CTT','CTC','CTA','CTG'],
'S' => ['TCT','TCC','TCA','TCG','AGT','AGC'],
'Y' => ['TAT','TAC'],
'C' => ['TGT','TGC'],
'W' => ['TGG'],
'P' => ['CCT','CCC','CCA','CCG'],
'H' => ['CAT','CAC'],
'Q' => ['CAA','CAG'],
'R' => ['CGT','CGC','CGA','CGG','AGA','AGG'],
'I' => ['ATT','ATC','ATA'],
'M' => ['ATG'],
'T' => ['ACT','ACC','ACA','ACG'],
'N' => ['AAT','AAC'],
'K' => ['AAA','AAG'],
'V' => ['GTT','GTC','GTA','GTG'],
'A' => ['GCT','GCC','GCA','GCG'],
'D' => ['GAT','GAC'], 
'E' => ['GAA','GAG'],
'G' => ['GGT','GGC','GGA','GGG'],
'X' => ['XXX'],
);

my @aa = ('A','R','N','D','C','Q','E','G','H','I','L','K','M','F','P','S','T','W','Y','V','X');

my $MAXPUNTUACION = 0; # máximo valor de la matriz de sustitución
                       # es una matriz de costes, con números negativos

my ($aa1,$aa2,$cod1,$cod2,@codones1,@codones2);
my ($dist,$distmedia);

print "#  Matriz de sustitucion de aminoacidos basada en el codigo genetico\n";
print "   A  R  N  D  C  Q  E  G  H  I  L  K  M  F  P  S  T  W  Y  V  X\n";
	
foreach $aa1 (@aa)
{
	@codones1 = @{ $tripletes{$aa1} };  # toma todos los codondes correspondientes a $aa1
	print "$aa1 ";
	
	foreach $aa2 (@aa)
	{
		if($aa1 eq $aa2){ print " $MAXPUNTUACION "; }
		else
		{
			@codones2 = @{ $tripletes{$aa2} };   # toma todos los codondes correspondientes a $aa2
			$distmedia = 0;
		
			foreach $cod1 (@codones1)
			{
				foreach $cod2 (@codones2)
				{
					# cuenta los cambios que hay que hacer para pasar de un codon a otro (max 3)
					$dist = 0;
					if(substr($cod1,0,1) ne substr($cod2,0,1)){ $dist++; }
					if(substr($cod1,1,1) ne substr($cod2,1,1)){ $dist++; }
					if(substr($cod1,2,1) ne substr($cod2,2,1)){ $dist++; }
					
					# recuerda la distancia entre estos dos codones para calcular la media entre aa1 y aa2 
					$distmedia += $dist;
				}
			}
			
			# calcula la distancia media entre aa1 y aa2
			$distmedia /= (scalar(@codones1)*scalar(@codones2));
			$distmedia = int($distmedia);
			
			# imprime el coste medio para pasar de aa1 a aa2
			print "-$distmedia ";
		}
	}
	print "\n";	
}

¿Se te ocurre cómo modificar el algoritmo?

Bruno Contreras-Moreira
http://www.eead.csic.es/compbio