Foro - Perl en Español

Hola perleros,

estoy interesado en crear un programa que me lea un archivo de secuencias cortas en formato FASTA y que me cuente el número de apariciones de cada secuencia. Me imagino que para ello tengo que crear un hash. La cuestión es que el archivo es muy grande y en pruebas con otro software he visto que puedo tener hasta 1 millón de secuencias únicas algunas repetidas más de una vez. Me gustaría luego poder crear un archivo de salida también en formato FASTA con los encabezados originales sólo para aquellas secuencias que aparecen al menos un número determinado de veces, especificado en los parámetros del programa. ¿Cómo puedo diseñar este programa? Me imagino que luego vendré con preguntas específicas, pero para empezar...

Primero, he pensado en crear un hash, meter la primera secuencia del archivo de entrada y luego empezar a comparar la siguiente. Si es igual, sumar uno al valor del key, si no, añadir otro key e inicializar su valor.

Después, he pensado eliminar todas aquellas entradas del hash que tengan un valor menor al especificado en los parámetros del programa.

La duda que me queda es cómo recuperar el encabezado original de cada uno de los keys que no se han eliminado y escribirlo a un archivo de salida.

Os incluyo un ejemplo de entrada:


Como veis, las secuencias 1 y 2 son iguales y las secuencias 3,4,5,6 y 7 también. Las tres últimas secuencias no se repiten.

La salida del programa, seleccionando un criterio de al menos 3 repeticiones debería ser:


¡Mil gracias!

D.

danusol escribiste:Primero, he pensado en crear un hash, meter la primera secuencia del archivo de entrada y luego empezar a comparar la siguiente. Si es igual, sumar uno al valor del key, si no, añadir otro key e inicializar su valor.

Después, he pensado eliminar todas aquellas entradas del hash que tengan un valor menor al especificado en los parámetros del programa.

Con los hash, no te hace falta "buscar" claves: es implícito en un hash.

La siguiente solución da la salida que pides. Observa que solo se trata de guardar las secuencias como las claves del hash, como indicabas, pero, en lugar de contar el número de veces que se repiten, guardamos, como valor de ellas, los nombres de las secuencias que las referencian. Tenemos entonces un hash que contiene arrays, y cada uno, es una lista de esos nombres.

Luego, a la hora de imprimir, deshacemos la madeja, quedándonos solo con los que cumplan la condición que buscamos.

La línea 24 es especial: este programa sólo funcionará si la secuencia está contenido en una sola línea, y que además sea la línea que sigue a la línea del nombre de la secuencia.

Explorer, muchísimas gracias por tu rápida respuesta. Ni se me había pasado por la cabeza la posibilidad de meter como valor un array, ¡increíble! (para un novato, claro). Es una solución muy inteligente. El problema que veo (y eso me pasa por no contar toda la historia) es que tengo en mente una modificación de este programa para que funcione con un formato similar, formato que en vez de dos líneas por secuencia, tiene cuatro, y ahora el indicador de nuevo bloque (cuatro líneas) es un "@", como sigue:



y entonces en este caso no se podría meter la primera, tercera y cuarta líneas dentro del valor de la clave en un array, ¿no?.

Por otro lado, ¿podrías explicarme un poco la expresión regular que usas para el encabezado? La primera parte la entiendo pero no me queda claro el significado de /x al final de la expresión en la línea 22.

Es normal, en Perl, pensar en "más de una dimensión" Pero solo hasta cierto punto, el suficiente para organizar la información que recibimos. Ideal, también, haber estudiado/leído algo de Estructuras de datos.

En lo referente al cambio de formato, pues entonces, hay que describirlo con detalle, cada una de sus partes. Por ejemplo, no sé la importancia de la línea que comienza por '+'. Es la secuencia la que hace de clave del hash, pero, como valor, puedes componer lo que quieras (la primera, la tercera y cuarta) en un solo escalar, que ese será el valor que meterás en el array (una posible solución).

Algo así:

Así, estamos metiendo las tres líneas en un solo valor escalar (con '.' concatenamos cadenas de caracteres).

El '/x' sirve para escribir la expresión regular con espacios en blanco, para que sea un poco más clara. Así, es lo mismo escribir

/^ > (.+)/x

que

/^>(.+)/

pero, un poco menos "apretado". Más información en perldoc perlre.

Muchísimas gracias de nuevo,

este es un formato de secuenciación masiva y las cuatro líneas son importantes: la primera es el encabezado del bloque; la segunda, la secuencia; la tercera, el encabezado de calidad de dicha secuencia; y la cuarta, la calidad de cada base en formato ASCII. La solución de concatenar las líneas 1, 3 y 4 me parece buena, pero tendré que añadir un guión bajo '_' o algo si quiero al final deconcatenar y poder recuperarlas y generar el archivo de resultados con el mismo formato de inicio.

Por cierto, efectivamente hablaba del valor de la clave para guardar el array y no de la clave en sí.

D.

No te hace falta unirlas con '_'... ya están separadas entre sí con el carácter "\n" que traen del mismo fichero del que fueron leídas.

Muchas gracias explorer, me pongo a ello. Aprendo muchísimo con el foro.

D.

Hola de nuevo.

Reactivo esta viejo hilo por estar muy relacionado con mi problema actual.

Tengo un archivo de texto con dos columnas separadas por espacio, de tal manera que la primera es un Transcript ID y la segunda es una proteína. El mismo Transcript ID puede corresponder a varias proteínas pero también puede haber el mismo Transcript ID con la misma proteína, por lo que existe redundancia, de la siguiente manera:



Lo que quiero hacer es ver cuántos Transcript ID representan cada proteína, pero sin redundancia, es decir si un mismo Transcript ID representa la misma proteína dos veces, entonces solo cuenta como uno.

Así, en el ejemplo de arriba, el resultado sería


Para ello, estoy intentando leer cada línea para formar un hash de array, de tal manera que cada clave sea una proteína y su valor sea el array de Transcript ID. Por último, quiero crear un array que me guarde los Transcript ID únicos en cada valor/array con la función 'uniq' y al final imprimir cada clave con la longitud de su correspondiente array único.

El resultado que obtengo no es el esperado, porque todos me da como resultado que tienen una única aparición.



¿Dónde fallo?

Gracias,

D.

¡ayyyy!

¡Mil gracias, ahora sí!

D.