|
0
|
1 #!/usr/bin/perl -w
|
|
|
2 $| = 1;
|
|
|
3
|
|
|
4 # Dado un fichero en formato GFF3 que incluye el análisis de varias secuencias,
|
|
|
5 # el programa devuelve un fichero de texto que incluye los valores comunes a todas las secuencias dadas
|
|
|
6 # para los atributos:
|
|
|
7 # Name -> Entrada de la base de datos de donde se ha obtenido una característica determinada;
|
|
|
8 # Ontology_term -> Entradas de Gene Ontology para una característica dada;
|
|
|
9 # Dbxref -> Entra de Interpro para una característica dada.
|
|
|
10
|
|
|
11 use strict;
|
|
|
12
|
|
|
13
|
|
|
14 # Declaración e inicialización de variables
|
|
|
15
|
|
|
16 my $fichero_ent = ""; #Nombre de fichero en formato GFF3 a analizar tomado de líneas de comandos
|
|
|
17 my $output=""; #Fichero de salida pasado como parámetro
|
|
|
18 my (@id_1, @id_2) = (); #Almacenan temporalmente las líneas directivas e ID de las secuencias
|
|
|
19 my @ids = (); #Todos los ID-seq del archivo
|
|
|
20 my @temp = (); #Todas las líneas de características del archivo
|
|
|
21 my $lin = ""; #Recupera cada ID del @ids
|
|
|
22 my @lin_id = (); #Todas las líneas de características para un ID determinado
|
|
|
23 my (@t1, @t3) = (); #Almacenan temporalmente las características y atributos de una línea dada
|
|
|
24 my $atributo = ""; #Únicamente la característica "Atributos" de cada línea de características
|
|
|
25 my @etiquetas = ("Name","Ontology_term","Dbxref"); #Son los tres tipos de atributos comunes que se van a extraer del fichero
|
|
|
26 my @sel_atrib = (); #Los atributos correspondientes a una etiqueta dada en cada fila
|
|
|
27 my ($etiq, $atrib) = ""; #La etiqueta y el valor respectivamente, de un atributo en una línea
|
|
|
28 my @val_atrib = (); #Los diferentes posibles valores de un atributo en una línea
|
|
|
29 my (@valores, @valores_rep) = ([],[],[]); #Valores de cada atributo para cada ID-seq y los valores repetidos para todos los ID-seq
|
|
|
30 my @repetidos = (); #Almacena temporalmente los valores comunes para cada atributo, entre los ID-seq analizados
|
|
|
31
|
|
|
32
|
|
|
33 ######## Abrir fichero y seleccionar lineas #########
|
|
|
34
|
|
|
35 ($fichero_ent,$output) = @ARGV;
|
|
|
36
|
|
|
37 open(ARCHGFF3, $fichero_ent) || die "Failure to open the file \"$fichero_ent\"\n\n"; #Abre el fichero
|
|
|
38
|
|
|
39 while (<ARCHGFF3>) { #Lee el archivo
|
|
|
40 chomp $_;
|
|
|
41 if ($_ =~ /^##FASTA/) { #Elimina la parte de secuencias fasta
|
|
|
42 last;
|
|
|
43 }elsif ($_ =~ /^##sequence-region./) {
|
|
|
44 push (@id_1, $_) #Las líneas directivas de sequence-region, para obtener su ID correspondiente
|
|
|
45 }elsif ($_ =~ /^#+/) { #Elimina líneas de comentarios y directivas, excepto el tipo anterior
|
|
|
46 next;
|
|
|
47 }else { push (@temp, $_)}; #Almacena las líneas con los atributos de todos los ID
|
|
|
48 };
|
|
|
49 close ARCHGFF3;
|
|
|
50
|
|
|
51 ########### Seleccionar ID de la línea directiva ##################
|
|
|
52
|
|
|
53 foreach (@id_1) {
|
|
|
54 @id_2 = split(/\s/,$_,3);
|
|
|
55 push (@ids, $id_2[1]); #Almacena todos los ID-seq que hay en el archivo
|
|
|
56 }
|
|
|
57
|
|
|
58 ########### Seleccionar un ID-seq determinado y todas sus filas de características ##################
|
|
|
59
|
|
|
60 my $i = 0; #Para distinguir entre la primera secuencia y el resto
|
|
|
61 foreach $lin (@ids) {
|
|
|
62 @lin_id = grep (/^$lin/, @temp); #Todas las líneas de características correspondientes al ID seleccionado
|
|
|
63
|
|
|
64 ########### Seleccionar, para cada línea de características, la columna novena de atributos ##################
|
|
|
65
|
|
|
66 foreach (@lin_id) {
|
|
|
67 @t1 = split(/\t/,$_); #Cada elemento es una característica de la línea dada
|
|
|
68 $atributo = $t1[8]; #Únicamente la característica "atributos" de la línea dada
|
|
|
69
|
|
|
70 #Seleccionar los atributos "Name", "Ontology_term" y "Dbxref" de la columna 9
|
|
|
71
|
|
|
72 @t3 = split(/;/, $atributo); #Cada elemento es un atributo de la característica "Atributos" de una línea dada
|
|
|
73
|
|
|
74 ########### Almacenar los diferentes valores de cada atributo (Name, Ontology_term y Dbxref) #############
|
|
|
75 ########### de un ID_seq determinado en un @rray diferente #############
|
|
|
76
|
|
|
77 for my $cont (0..2){
|
|
|
78 if (@sel_atrib = grep (/^$etiquetas[$cont]./, @t3)) { #Evitar valores no definidos para un atributo concreto
|
|
|
79 ($etiq, $atrib) = split (/=/, $sel_atrib[0], 2);
|
|
|
80 @val_atrib = split (/,/, $atrib);
|
|
|
81 foreach my $valor (@val_atrib) {
|
|
|
82 if (!grep (/$valor/, @{$valores[$cont]})) { #Evitar valores repetidos de un mismo atributo para un ID-seq determinado
|
|
|
83 push (@{$valores[$cont]}, $valor); #Todos los valores diferentes de cada atributo para un ID-seq determinado
|
|
|
84 }
|
|
|
85 }
|
|
|
86 }
|
|
|
87 }
|
|
|
88 }
|
|
|
89
|
|
|
90 ########### Comprobar los atributos comunes a todos los ID-seq y guardarlos #################
|
|
|
91
|
|
|
92 if ($i == 0) { #Para el primer ID-seq se guardan todos sus atributos
|
|
|
93 for my $cont (0..2) { #no repetidos. Como mucho, serán todos estos la solución.
|
|
|
94 $valores_rep[$cont] = [@{$valores[$cont]}];
|
|
|
95 $valores[$cont] = [];
|
|
|
96 }
|
|
|
97 $i++;
|
|
|
98 } else {
|
|
|
99 for my $cont (0..2) {
|
|
|
100 foreach my $valor (@{$valores[$cont]}) {
|
|
|
101 if (grep (/$valor/, @{$valores_rep[$cont]})) { #Búsqueda valores comunes para cada tipo atributo entre los
|
|
|
102 push (@repetidos, $valor); #ID-seq analizados.
|
|
|
103 }
|
|
|
104 }
|
|
|
105 $valores_rep[$cont] = [@repetidos]; #Todos los valores comunes de los tres tipos de atributos para
|
|
|
106 $valores[$cont] = []; #todas las secuencias del fichero.
|
|
|
107 @repetidos = ();
|
|
|
108 }
|
|
|
109 }
|
|
|
110
|
|
|
111 }
|
|
|
112
|
|
|
113
|
|
|
114 ############ Impresión de resultados ######################
|
|
|
115
|
|
|
116 if (!open(FICHEROUT, ">$output")) {
|
|
|
117 print "The file \"$output\" can not be opened";
|
|
|
118 } else {
|
|
|
119 print FICHEROUT "The common atributes of the sequences ","@ids"," are:\n";
|
|
|
120 for my $cont (0..2) {
|
|
|
121 print FICHEROUT "For atribute ","$etiquetas[$cont]" ," : ";
|
|
|
122 if (scalar @{$valores_rep[$cont]} == 0) {
|
|
|
123 print FICHEROUT "No common atributes \n";
|
|
|
124 }else { print FICHEROUT "@{$valores_rep[$cont]}", "\n"}
|
|
|
125 }
|
|
|
126 close (FICHEROUT);
|
|
|
127 }
|
|
|
128 exit;
|
|
|
129
|
|
|
130
|
|
|
131
|
|
|
132
|
