Alineación de múltiples secuencias en un grupo heterogéneo

Tengo una lista de ~200 secuencias de ADN, que representan probablemente 50 regiones genómicas diferentes, pero todas están mezcladas. Por ejemplo, si tengo seq1, seq2.... seq10, seq1podría alinearse con seq3y seq8, pero no tener ninguna relación con los demás.

También hay algunas diferencias en las longitudes muestreadas, por lo que el ejemplo anterior podría representar:

Seq1-------------------------------------------------

Seq3---------------------- seq8-----------------

Tal que seq3y seq8no se alinean entre sí, pero ambos se alinean paraseq1

Entonces, lo que me gustaría hacer es revisar y de alguna manera generar una lista de grupos de secuencias que se alinean entre sí, así como las alineaciones. p.ej:

Grupo 1

Seq1-------------------------------------------------

Seq3---------------------- seq8-----------------

Grupo 2

seq2------------------------------------

. seq6-----------------------------

. seq7--------------------------xxxxxx

Grupo 3 ... etc.

Intentar ClustalW o MUSSLE para alinear todo no funciona (o lleva una cantidad de tiempo irrazonable), supongo que porque hay muchas secuencias que no se alinean en absoluto. Intenté crear una base de datos BLAST personalizada y luego BLASTé cada secuencia contra ella, pero luego obtengo múltiples aciertos para la misma alineación (con el ejemplo del grupo 2 anterior, , , , , seq2:seq6y seq2:seq7todos seq6:seq2se seq6:seq7devuelven seq7:seq2como seq7:seq66 aciertos únicos, cuando deberían ser agrupados.

Mi conocimiento de codificación actual es bastante básico, pero estoy dispuesto a leer documentos y resolver cosas, simplemente no quiero reinventar la rueda.

Edit2: Realmente, la agrupación es la parte importante: una vez que tengo los grupos, puedo hacer la alineación por separado con poco esfuerzo. Solo me gustaría tener grupos donde cada secuencia esté en un solo grupo.

Respuestas (2)

Puedes probar esto:

  • BLAST cada secuencia a cada otra secuencia (por pares).
  • Cada alineación (con algún corte definido) denota una conexión.
  • Mapear todas las conexiones.
  • Si una secuencia está conectada con otra directa o indirectamente, cae en un grupo. Coloque todas las secuencias que seq1se alinean en el Grupo-1 , luego vaya a las alineaciones de estas secuencias; poner todas las secuencias a las que se alinean de nuevo en el Grupo-1 ; así que sigue poblando el grupo de esta manera.

Metodología:

  • Instale blast independiente (si no tiene muchas secuencias, también puede ejecutar BLAST en línea)
  • Cree una base de datos explosiva a partir de sus secuencias usandomakeblastdb
  • Alinee estas secuencias con la base de datos. Si está utilizando BLAST en línea, use BL2seq (alinee dos secuencias). Es mucho mejor y conveniente usarlo de forma independiente. También puede mencionar si desea más-más o más-menos o ambas alineaciones. En algunos casos, es posible que desee solo cualquiera de los dos.
  • En el BLAST independiente, puede especificar el formato de salida (qué campos incluir, etc. El formato que elija solo depende de sus requisitos).

Un formato de salida tabular se parece a esto:

# BLASTN 2.2.27+
# Query: TCONS_00036712 gene=XLOC_017996
# Database: ../nt_db/nt
# Fields: query id, subject id, % identity, alignment length, mismatches, gap opens, q. start, q. end, s. start, s. end, evalue, bit score
# 1014 hits found
TCONS_00036712  gi|191174875|emb|CU655970.6|    95.54   202 9   0   423 624 16680   16479   8e-85    324
TCONS_00036712  gi|51491599|gb|AC144709.2|  95.02   201 10  0   424 624 28443   28243   1e-82    316

Ignore las #líneas comentadas ( ); el primer campo es el ID de la consulta, el segundo es el ID del asunto y hay una alineación entre los dos; los otros campos brindan información sobre la alineación (puede elegir estos campos).

Para el análisis utilizo un lenguaje de secuencias de comandos rápido y fácil llamado awk que se incluye en todos los sistemas basados ​​en UNIX. También está disponible para Windows en el paquete GNUWin32.

Lo que tienes que hacer es revisar los dos primeros campos y actualizar el grupo.

# CrearGrupos.awk

BEGIN{FS="\t"} # Declarando el separador de campo como tabulador

!($1 en grp){ # Comprobar si la secuencia es un grupo principal. Si no...
    k=1
    para(yo en grp){
        if($1 in grp[i]){ # Comprobar si la secuencia es parte de algún otro grupo
            parentgrp[$1]=i
            if(!($2 in grp[i])) # Comprobar si el segundo campo, es decir, el sujeto, ya está presente en el grupo principal
                grp[i][$2] # si no asigna el segundo campo al grupo principal
            k=0
            romper # dejar de comprobar más
        }
    }
    if(k==1) # No hay un grupo padre con esta etiqueta y la secuencia no es parte de ningún otro grupo
        grp[$1][$1] # Cree un grupo con el id de consulta como etiqueta y agregue esta consulta a este grupo.
}  

$1 en grp{
    if(!($2 en grp[$1]))
        grp[$1][$2] # Si el segundo campo no es parte del grupo con el primer campo como etiqueta, asígnelo a ese
}

FIN{
    para(yo en grp){
        x++
        imprimir "Grupo-"x"\n----------"
        for(j en grp[i])
            imprimir j
    }
    imprimir "\n"
}

Ejecute este script como este en la terminal:
gawk -f MakeGroups.awk blastalignmentfile.txt

Nota: este script contiene matrices multidimensionales. No funcionará con todas las versiones de awk. usogawk version >4.0 _

Como menciona swarnbes en su respuesta, existen algoritmos más rápidos que hacen este tipo de cosas y se utilizan para el ensamblaje de secuencias. Lo que hacen muchos de ellos es hacer un gráfico (redes llamadas gráficos de deBruijn), donde cada conexión es una alineación, y calcular un Camino Euleriano. Vea esta reseña de Pavel Pevzner para más detalles. Las secuencias superpuestas forman contigs y puede rastrear fácilmente qué secuencia proviene de qué contig (que puede llamar Grupo). Cada contig/Group es un subgrafo deBruijn disjunto.

Sí, aunque si pudieras señalarme la dirección correcta, sería genial.
@kevbonham Este script debería funcionar. Si te gusta algún otro lenguaje de programación, puedes implementar esta lógica en eso.
Wow, muchas gracias! Parece que debería funcionar... Podría intentar volver a escribir en python (ya que es con lo que estoy más familiarizado), pero la lógica es bastante fácil de seguir. Supongo que agregaría otro bucle para ignorar las coincidencias <un cierto % de coincidencia, ¿o especificaría esos parámetros en los parámetros iniciales de BLAST?
@kevbonham puede agregar esa condición junto con verificar si el segundo campo está en un grupo o no, sin necesidad de agregar otro ciclo. Básicamente, todos esos criterios deben ser estrictos con la relación entre el primer y el segundo campo. Lo que pasa con awk es que transmite lecturas (no carga todo el archivo en la RAM) y la separación de campos es un trabajo fácil. Bueno para el análisis de texto.
Tiene sentido. ¿La lectura de secuencias es equivalente a algo como readline()declarado explícitamente en python?
@kevbonham readline()lee el archivo completo línea por línea y lo almacena como una lista. En la lectura continua, se lee una línea y se implementa el código en ella antes de pasar a la línea siguiente. Para implementar la lectura de secuencias en el uso de python for line in file:. Ver aquí _ Awk por defecto hace lectura de flujo. Awk también tiene una buena base de expresiones regulares; en python tienes que hacerlo import re. Dicho esto, Awk no es tan versátil como Python, pero es bueno para ciertas tareas.
@kevbonham en awk puede especificar cuál es el separador de registros (RS, que es una nueva línea de forma predeterminada). En su lugar, puede establecer ">.*"como RS para un archivo fasta y cada secuencia de nucleótidos/proteínas se leerá en una secuencia.
@kevbonham Otro punto: el equivalente de matrices asociativas de awk (hashtables) en python es un diccionario. Los diccionarios a veces son muy convenientes y funcionan mucho más rápido que las listas.
Finalmente logré hacer que esto funcione en python, usando su script como plantilla (al menos lo que entendí de la lógica detrás de esto). Probablemente esté mucho más hinchado de lo necesario, pero fue un buen ejercicio de aprendizaje. ¿Es apropiado publicar mi código como respuesta a esta pregunta como referencia? Todavía no estoy completamente seguro acerca de la etiqueta de intercambio de pilas.
@kevbonham Sí, puede hacerlo, pero esto no debería tratarse de escribir códigos en general y este no es el foro adecuado para profundizar en los detalles de la programación. El código que explica un algoritmo está bien. El pseudocódigo es mejor. La escritura de códigos se encuentra bajo el dominio de StackOverflow.

¿Realmente necesitas BLAST? Es decir, ¿las secuencias son lo suficientemente diferentes entre sí como para necesitar un algoritmo que busque grandes diferencias entre ellas?

Tal vez podría usar algo como Phrap, que debería armar contigs para usted, si las secuencias que deberían ir juntas son muy parecidas.

Bueno no exactamente. Donde las secuencias se superponen, se superponen con >99 % de identidad, pero puede haber secuencias laterales que no se superponen. Cambié un poco la pregunta para aclarar (espero).
No estoy seguro de si Phrap en sí mismo funcionaría porque está diseñado para el ensamblaje de secuencias, pero el cross_match o SWAT relacionado podría funcionar... los investigaré un poco más.
Alineación de múltiples secuencias en un grupo heterogéneo
RespuestasAquíPolítica de privacidad1y, 0v