Foro - Perl en Español

Bueno, esta es mi solución, en Perl moderno.

Buenos días pvaldes y explorer. Muchas gracias por vuestra ayuda y sugerencias. Vamos por partes.

@pvaldes, los genes duplicados son una desgracia pero los necesito. Esto es Arabidopsis y existen esos genes en distintos cromosomas con su ID único pero con mismo gene_symbol. Por lo tanto a un mismo gene_symbol le pueden corresponder uno o más ID. Pero a un ID le corresponde un único gene_symbol. Los gene_symbol se pueden duplicar pero los ID no. Por eso la necesidad de este script, para que cuando tengo mi lista con gene_symbols y posición en el genoma saber a qué ID se refiere porque conozco su posición. Mi resultado final tiene que ser donde tiene que haber tantos gene_symbol como había en la "lista" inicial, y con su correspondiente ID, para luego poder seguir trabajando con el ID.

@explorer, qué puedo decir, eso es música. De momento como ya he dicho anteriormente, suprimiría la línea 28 porque no necesito saber qué hay en cada cromosoma, Eso ya viene de partida en la lista y en la DB. Por otro lado, me sorprende que el "split" sea tan inteligente que no le tenga que decir cuál es el separador. De hecho quiero que me separe únicamente por tabulación porque algunos nombres tienen algún punto o espacio y no quiero que me los separe. Estoy intentando decirle que me separe por tabulación y no me lo hace:



no sería así?

Por otro lado, la estructura creada para almacenar la información de mi lista la tengo que estudiar porque eso de un hash de arrays de arrays ¡me crea cortocircuito! Por último el bucle for() es muy curioso. También tengo que estudiármelo. Voy a probar que esto funcionaría con mi lista total y mi base de datos total.



Por cierto, ¿qué ventaja o diferencia hay entre usar "say" y "print"?

Un posible problema podría ser el caso que hubiese dos genes en la base de datos con la misma posición cromosómica pero por ejemplo uno en lectura "sentido" y otro en lectura "antisentido". Serían dos genes con distinto $symbol, distinto $ID pero similar $start y $end. Entonces mi $pos cumpliría los requisitos para ambos y el programa me sacaría dos líneas de resultado. ¿Cómo podría meter en ese último 'for' que una vez que encuentra un $pos que cumple la condición, pase al siguiente $pos y no siga buscando el mismo en el resto de la base de datos? Esto podría ser un problema, lo sé, habría que incluir en el programa información de si está en "sentido" o "antisentido", pero ¿y si me diese igual y cojo el primer "acierto" como válido?

Un saludo,

D.

danusol escribiste:De momento como ya he dicho anteriormente, suprimiría la línea 28 porque no necesito saber qué hay en cada cromosoma.

Es obvio que las líneas 3 y 28 solo se usan para mostrar el aspecto de la estructura de datos tridimensional que vamos a usar. Es solo un visionado. Cuando se termina de entender, se comentan o eliminan esas líneas.

danusol escribiste:me sorprende que el "split" sea tan inteligente que no le tenga que decir cuál es el separador.

split() sin argumentos, parte a la variable $_ usando como delimitador cualquier conjunto de espacios en blanco. Más información en perldoc -f split.

danusol escribiste:la estructura creada para almacenar la información de mi lista la tengo que estudiar porque eso de un hash de arrays de arrays ¡me crea cortocircuito!

¡Bah!, ¡no tanto! Solo son tres dimensiones... Mi récord son 8 dimensiones

danusol escribiste:Por último el bucle for() es muy curioso. También tengo que estudiármelo.

En el bucle recorremos todas las posiciones para el mismo cromosoma. Por cada $posición, analizamos:

* Si la $posición es menor que el $inicio, hacemos un next para que analice la siguiente posición, pues eso significa que aún no hemos entrado en el rango [$inicio,$final] (en notación matemática)

* Si la $posición es mayor que el $final, es que nos hemos pasado del rango, así que terminamos inmediatamente con un last. No es necesario seguir mirando más posiciones, pues, al estar ordenadas (líneas 24 a 26), ya sabemos que las restantes seguirán quedando fuera del rango.

* Si llega a la línea 50, es que se trata de una $posición que está en el rango [$inicio,$final], así que lo imprimimos.

danusol escribiste:Por cierto, ¿qué ventaja o diferencia hay entre usar "say" y "print"?

Nos ahorramos poner un "\n" al final. Más(?) información en perldoc -f say.

explorer escribiste:¡Bah!, ¡no tanto! Solo son tres dimensiones... Mi récord son 8 dimensiones

entonces seguro que entendiste "Origen" mucho mejor que yo, que me costó situarme al final!

explorer escribiste:En el bucle recorremos todas las posiciones para el mismo cromosoma. Por cada $posición, analizamos:

* Si la $posición es menor que el $inicio, hacemos un next para que analice la siguiente posición, pues eso significa que aún no hemos entrado en el rango [$inicio,$final] (en notación matemática)

* Si la $posición es mayor que el $final, es que nos hemos pasado del rango, así que terminamos inmediatamente con un last. No es necesario seguir mirando más posiciones, pues, al estar ordenadas (líneas 24 a 26), ya sabemos que las restantes seguirán quedando fuera del rango.

* Si llega a la línea 50, es que se trata de una $posición que está en el rango [$inicio,$final], así que lo imprimimos.

Desgraciadamente estoy encontrando $posicion que está en el rango [$inicio,$final] de dos IDs. Entonces al ejecutar el programa el resultado me da que $gen sale dos veces. De hecho si mi "lista" inicial tiene 2703 entradas, mi resultado tiene 2732 líneas. He encontrado el siguiente ejemplo:



por lo que un mismo $pos es seguro que valida la condición de que estén entre estos dos rangos. Parece que puede darse el caso, posiblemente por la dirección del gen. Necesito que sólo me dé un resultado, aunque puede ser que no sea el cierto (luego al estudiar los resultados podría dar con el real) y que mi archivo de salida tenga las mismas entradas que mi "lista".

En el ejemplo que pusiste, el gen ACO2 sale dos veces, pero está en posiciones distintas, por lo que es correcto que salga, pero lo que quieres es que no salga si coinciden las posiciones, ¿no?

Y según las reglas que has puesto, la salida es posible que no tenga las mismas líneas que la "lista": los genes con posiciones que están fuera de los rangos, no salen.

¿No sería así?

Depende de los tabuladores que tengas, puede encontrar problemas si hay más de dos.

Una opción posible si sabes que el problema está en el último elemento es darle más libertad y dejar que se encargue de recoger todo


Aquí, @gen funcionaría como un resto de número variable de elementos. Un join() posterior permitiría devolverlo a una cadena con rapidez así que directamente capturamos todos los nombres.

Otra opción posible es añadir un límite al número de campos posibles para trocear, en éste caso 3

pvaldes escribiste:

Sintáxis: [ Descargar ] [ Ocultar ]

Sintáxis: [ Descargar ] [ Mostrar ]

Using perl Syntax Highlighting

1. my ($cromosoma, $posición, $gen) = split /\t/;

Coloreado en 0.001 segundos, usando GeSHi 1.0.8.4

¿No sería así?

Depende de los tabuladores que tengas, puede encontrar problemas si hay más de dos.

Estrictamente, y viendo el código original sí, pero...

split(), sin argumentos, corta $_ por los espacios en blanco o tabuladores, sea el número que sea, uno, dos, cien...

Por eso, aparte de quedar bonito en el código, solo poner split, también nos vale por si la entrada cambia (los tabuladores, un día, cambian a 8 espacios)

explorer escribiste:En el ejemplo que pusiste, el gen ACO2 sale dos veces, pero está en posiciones distintas, por lo que es correcto que salga, pero lo que quieres es que no salga si coinciden las posiciones, ¿no?

Y según las reglas que has puesto, la salida es posible que no tenga las mismas líneas que la "lista": los genes con posiciones que están fuera de los rangos, no salen.

Explorer, disculpa por el lío. A ver si consigo explicarlo mejor.

Partimos de que todos los genes de "lista" deben estar en "DB". "lista" es el resultado de un experimento en Arabidopsis y "DB" es la base de datos completa de genes de Arabidopsis.

Ahora bien, en "lista" hay genes, como ACO2, pero soy capaz de identificarlos automáticamente en "DB" porque tienen distinto cromosoma y distinta posición. Por eso ACO2 no me da problemas porque están en distinta posición.
Pero en "DB" hay IDs distintos cuyo [$inicio,$final] es prácticamente igual y el cromosoma es el mismo. Puede darse el caso por tanto que una posición de "lista" cumpla las condiciones de dos ID de "DB". Te intento poner ejemplos de los dos casos:

En lista podemos tener esto

Chr5 25894017 ACL5
Chr2 8476733 ACO1
Chr4 13546819 ACO2
Chr1 23083211 ACO2
Chr2 1138898 ACR5

En DB podemos tener esto

AT5G64750 Chr5 25891505 25896528
AT5G64760 Chr5 25891505 25896582
AT3G18710 6434083 6435561 Chr3
AT4G25880 13154708 13159382 Chr4
AT1G71695 26964249 26966688 Chr1
AT5G41480 16595927 16598636 Chr5
AT5G15008 4856975 4857139 Chr5
AT3G18310 6284356 6287144 Chr3
AT1G62380 23082170 23084253 Chr1
AT4G26970 13542918 13548629 Chr4

Para ACO2 en azul, el programa lo asignará AT4G26970
Para ACO2 en rojo, el programa lo asignará AT1G62380

pero en el caso de ACL5 en negrita, encuentra dos posibles resultados: AT5G64750 y AT5G64760 porque posición está dentro de los dos intervalos. Como no puedo averiguar cuál es automáticamente, prefiero asignarle el primero que aparezca y luego al estudiar con detalle todo veré que no era el correcto, que no asignarle los dos y acabar con unos resultados más grandes que la lista inicial. Por eso decía que si una posición de "lista", digamos posición1, cumple las condiciones en "DB" una vez, que pase a posición2 y no que siga buscando posibles resultados para posición1.

Entonces mi resultado tendría que ser

ACL5 AT5G64750
ACO1 AT_loquesea
ACO2 AT4G26970
ACO2 AT1G62380
ACR5 AT_loquesea



Espero haber clarificado el problema.

¡Anda!

¿Eso es correcto? ¿Pueden aparecer registros en que están cambiados de sitio el cromosoma y las posiciones?

No, lo siento, me he confundido al coger los datos ya que antes tengo que preparar la base de datos y la "lista" para reducirlas a la mínima información necesaria. Perdón por el error.

Los campos deben ser como hemos estado trabajando hasta ahora.
Para DB:
Por lo tanto sería correcto así:


Para la "lista" sería: