#!/usr/bin/perl -w
use strict;
#####################################################
#################### ####################
#################### INFORMACIO ####################
#################### ####################
#####################################################
#
### NOM DEL PROGRAMA: EXON-finding
#
### AUTORES: Alba Llop, Marta Mancebo i Mireya Plass
#
### DATA: 18 -03 - 2004
#
### 4t BIOLOGIA HUMANA, Universitat Pompeu Fabra
#
#####################################################
#####################################################
#################### ####################
#################### PRINCIPAL ####################
#################### ####################
#####################################################
##### DECLARACIO DE VARIABLES GLOBALS #############################################################################################
my $id; # emmagatzema el nom de la nostra sequencia.
my $sequencia; # guarda la sequencia en format fasta en un string.
my %donors; # emmagatzema la matriu dels donors en un hash.
my $pos1; # emmagatzema el numero de posicions que s'examinen amb la matriu de donors.
my %acceptors; # emmagatzema la matriu dels acceptors en un hash.
my $pos2; # emmagatzema el numero de posicions que s'examinen amb la matriu d'acceptors.
my %codons; # emmagatzema la matriu dels codons en un hash.
my $llindar; # emmagatzema el llindar que utilitzarem per determinar si un senyal (acceptor o donor) es bo o no.
my $global; # emmagatzema la puntuacio minima que volem que tingui un exo predit.
my $prediccio; # guarda el nombre d'exons predits que es vol obtenir per cada senyal (emparellaments).
my $tamany; # guarda la longitud minima que volem que tingui un exo predit.
my %ORF; # emmagatzema per a cadascuna de les pautes de lectura els diferents ORFs trobats.
my @exons; # guarda un array d'arrays on per cada exo guardem la posicio d'inici, la de final i la puntuacio i el frame.
my @first; # guarda els resultats amb els diferents exons First predits.
my @internal; # guarda els resultats amb els diferents exons Internal predits.
my @terminal; # guarda els resultats amb els diferents exons Terminal predits.
my @final; # guarda els resultats finals amb els diferents gens predits (inici, final, puntuació i frame).
my @resultats; # variable de treball. Guarda resultats de funcions que son arrays.
my $resultats; # variable de treball. Guarda resultats de funcions que son scalars.
###################################################################################################################################
##### OBRIM ELS ARGUMENTS QUE PASSEM A LA FUNCIO I ELS GUARDEM EN DIFERENTS VARIABLES #############################################
if (scalar(@ARGV) < 8) {
print "entrada de dades incorrecta.\nPosar fitxer sequencia fasta, fitxer matriu donors,\nfitxer matriu acceptors, fitxer matriu codons,\n llindar, nombre de prediccions,\n puntuacio global minima d'exo longitud minima d'exo\n";
exit(1);
}
# assignem en diferents variables els diferents arguments introduits per l'execucio del programa.
my $seq = $ARGV[0]; # guarda el nom del fitxer de la sequencia que hem d'obrir.
my $donor = $ARGV[1]; # guarda el nom del fitxer de la matriu de donors.
my $accep = $ARGV[2]; # guarda el nom del fitxer de la matriu dels acceptors.
my $codon = $ARGV[3]; # guarda el nom del fitxer de la matriu dels codons.
$llindar = $ARGV[4];
$prediccio = $ARGV[5];
$global = $ARGV[6];
$tamany = $ARGV[7];
###################################################################################################################################
##### CRIDEM LA FUNCIO form #######################################################################################################
@resultats = &form; # assignem els resultats de la funcio form a la variable @resultats.
$id = ${$resultats[0]}; # desreferenciem l'identificador i l'emmagatzemem en una variable singular.
$sequencia = ${$resultats[1]}; # desreferenciem la sequencia i l'emmagatzemem en un string.
###################################################################################################################################
##### CRIDEM LA FUNCIO mat ########################################################################################################
@resultats = mat($donor); # assignem els resultats de la funcio mat sobre l'arxiu que conte
# la matriu de pesos dels donors a la variable @resultats.
$pos1 = ${$resultats[0]}; # desreferenciem i assignem l'element que hi ha a la 1a posicio de @resultats.
%donors = %{$resultats[1]}; # desreferenciem i assignem l'element que hi ha a la 2a posicio de @resultats.
@resultats = mat($accep); # assignem els resultats de la funcio mat sobre l'arxiu que conte
# la matriu de pesos dels acceptors a la variable @resultats
$pos2 = ${$resultats[0]}; # desreferenciem i assignem l'element que hi ha a la 1a posicio de @resultats.
%acceptors = %{$resultats[1]}; # desreferenciem i assignem l'element que hi ha a la 2a posicio de @resultats
###################################################################################################################################
##### CRIDEM LA FUNCIO ORF ########################################################################################################
# cridem la funcio ORF i guardem el hash que conte els diferents ORFs trobats
# en una de les claus que identifica el frame en que el trobem.
$resultats = ORF($sequencia,0); # assignem la referencia dels resultats de la funcio ORF.
%{$ORF{0}} = %{$resultats}; # desreferenciem i asignem els resultats sl hash que hi ha dins de $ORF{0}.
$resultats = ORF($sequencia,1); # assignem la referencia dels resultats de la funcio ORF.
%{$ORF{1}} = %{$resultats}; # desreferenciem i asignem els resultats sl hash que hi ha dins de $ORF{1}.
$resultats = ORF($sequencia,2); # assignem la referencia dels resultats de la funcio ORF.
%{$ORF{2}} = %{$resultats}; # desreferenciem i asignem els resultats sl hash que hi ha dins de $ORF{2}.
###################################################################################################################################
##### CRIDEM LA FUNCIO pes ########################################################################################################
# guardem per cada frame (key 0, 1 o 2) un hash que indexa els
# diferents ORFs i que conte on array on s'emmagatzemen l'inici,
# el final, un hash amb els donors i un hash amb els acceptors.
$resultats = pes($sequencia, \%{$ORF{0}}, \%donors, $pos1, $llindar, 2); # assignem la referencia als resultats de la funcio pes.
%{$ORF{0}} = %{$resultats}; # desreferciem i guardem els resultats que hi ha a $ORF{0}.
$resultats = pes($sequencia, \%{$ORF{1}}, \%donors, $pos1, $llindar, 2); # assignem la referencia als resultats de la funcio pes.
%{$ORF{1}} = %{$resultats}; # desreferciem i guardem els resultats que hi ha a $ORF{1}.
$resultats = pes($sequencia, \%{$ORF{2}}, \%donors, $pos1, $llindar, 2); # assignem la referencia als resultats de la funcio pes.
%{$ORF{2}} = %{$resultats}; # desreferciem i guardem els resultats que hi ha a $ORF{2}.
$resultats = pes($sequencia, \%{$ORF{0}}, \%acceptors, $pos2, $llindar, 3); # assignem la referencia als resultats de la funcio pes.
%{$ORF{0}} = %{$resultats}; # desreferciem i guardem els resultats que hi ha a $ORF{0}.
$resultats = pes($sequencia, \%{$ORF{1}}, \%acceptors, $pos2, $llindar, 3); # assignem la referencia als resultats de la funcio pes.
%{$ORF{1}} = %{$resultats}; # desreferciem i guardem els resultats que hi ha a $ORF{1}.
$resultats = pes($sequencia, \%{$ORF{2}}, \%acceptors, $pos2, $llindar, 3); # assignem la referencia als resultats de la funcio pes.
%{$ORF{2}} = %{$resultats}; # desreferciem i guardem els resultats que hi ha a $ORF{2}.
##################################################################################################################################
##### CRIDEM LA FUNCIO def #######################################################################################################
@resultats = def(\%{$ORF{0}},$prediccio, $sequencia, 0, $tamany); # assignem les referencies als resultats de la funcio def.
@{$exons[0][0]} = @{$resultats[0]}; # guardem els arrays dels exons First del frame 0.
@{$exons[0][1]} = @{$resultats[1]}; # guardem els arrays dels exons Internal del frame 0.
@{$exons[0][2]} = @{$resultats[2]}; # guardem els arrays dels exons Terminal del frame 0.
@resultats = def(\%{$ORF{1}},$prediccio, $sequencia, 1, $tamany); # assignem les referencies als resultats de la funcio def.
@{$exons[1][0]} = @{$resultats[0]}; # guardem els arrays dels exons First del frame 1.
@{$exons[1][1]} = @{$resultats[1]}; # guardem els arrays dels exons Internal del frame 1.
@{$exons[1][2]} = @{$resultats[2]}; # guardem els arrays dels exons Terminal del frame 1.
@resultats = def(\%{$ORF{2}},$prediccio, $sequencia, 2, $tamany); # assignem les referencia als resultats de la funcio def.
@{$exons[2][0]} = @{$resultats[0]}; # guardem els arrays dels exons First del frame 2.
@{$exons[2][1]} = @{$resultats[1]}; # guardem els arrays dels exons Internal del frame 2.
@{$exons[2][2]} = @{$resultats[2]}; # guardem els arrays dels exons Terminal del frame 2.
###################################################################################################################################
##### CRIDEM LA FUNCIO bcod #######################################################################################################
$resultats = bcod($codon); # guardem la referencia als resultats de la funcio bcod.
%codons = %{$resultats}; # desreferenciem i guardem el hash dels codons a %codons.
###################################################################################################################################
##### CRIDEM LA FUNCIO biaix ######################################################################################################
@resultats = @{biaix($sequencia, \%codons, \@{$exons[0]}, $global)}; # guardem les referencies dels resultats de la funcio biaix
$first[@first] = $resultats[0]; # guardem els arrays dels exons First del exons del frame 0.
$internal[@internal] = $resultats[1]; # guardem els arrays dels exons Internal del frame 0.
$terminal[@terminal] = $resultats[2]; # guardem els arrays dels exons Terminal del frame 0;
@resultats = @{biaix($sequencia, \%codons, \@{$exons[1]}, $global)}; # guardem les referencies dels resultats de la funcio biaix
$first[@first] = $resultats[0]; # guardem els arrays dels exons First del frame 1.
$internal[@internal] = $resultats[1]; # guardem els arrays dels exons Internal del frame 1
$terminal[@terminal] = $resultats[2]; # guardem els arrays dels exons Terminal del frame 1.
@resultats = @{biaix($sequencia, \%codons, \@{$exons[2]}, $global)}; # guardem les referencies dels resultats de la funcio biaix
$first[@first] = $resultats[0]; # guardem els arrays dels exons First del frame 2.
$internal[@internal] = $resultats[1]; # guardem els arrays dels exons Internal del frame 2.
$terminal[@terminal] = $resultats[2]; # guardem els arrays dels exons Terminal del frame 2.
###################################################################################################################################
##### CRIDEM LA FUNCIO junta ######################################################################################################
@{$final[0]} = @{junta(\@first)}; # assignem l'array de resultats de la funcio junta a la posicio 0 de @final.
@{$final[1]} = @{junta(\@internal)}; # assignem l'array de resultats de la funcio junta a la posicio 1 de @final.
@{$final[2]} = @{junta(\@terminal)}; # assignem l'array de resultats de la funcio junta a la posicio 2 de @final.
@final = (@{$final[0]},@{$final[1]},@{$final[2]}); # guarda totes les prediccions juntes al mateix vector
###################################################################################################################################
##### ORDENEM ELS ARRAYS QUE HI HA AL HASH %final #################################################################################
@final = sort{${$a}[0] <=> ${$b}[0] || ${$a}[1] <=> ${$b}[1]}(@final); # ordenem els exons per la posicio d'inici
# i la posicio de final (1a i 2a columna).
###################################################################################################################################
##### FEM EL PRINT DE LES DADES EN FORMAT gff #####################################################################################
$id =~ s/>//; # eliminem el simbol ">" de l'identificador de la sequencia.
##### RESSEGUIM L'ARRAY QUE GUARDA ELS EXONS #####
my $c = 0; # inicialitzem el comptador de la iteracio a zero.
while ($c < @final){ # mentre hi hagi prediccions d'exons.
my @print = @{$final[$c]}; # assignem la informacio de l'exo al vector @print a cada volta.
print $id, "\t", "EXON-finding\t", $print[4],"\t", $print[0], "\t", $print[1], "\t";
printf "%.3f",$print[2];
print "\t","+\t",$print [3],"\n";
$c++; # incrementem el comptador de la iteracio.
}
###################################################################################################################################
###################################################################################################################################
####################################################
#################### ####################
#################### FUNCIONS ####################
#################### ####################
####################################################
# FUNCIO 1(form): obre l'argument que conte tot l'arxiu i retorna la referencia de l'identificador i
# la referencia de la sequencia directa en majuscules guardada en una string
#
# Parametres requerits: cap (la funcio obre directament el filehandle).
sub form {
my $s = ""; # guarda la sequencia en una variable singular (tota seguida) i la inicialitzem amb la paraula buida
my $id; # guarda l'identificador de la sequencia
my @d; # vector que guarda la seq directa
my $c = 0; # comptador de la iteracio
if (!open(SEQ,"< $seq")) { # obrim la sequencia emmagatzemada a l'argument
print " NO S'HA POGUT OBRIR EL FITXER DE LA SEQUENCIA!!!\n"; # i, en cas que no es pugui obrir,fem saltar el
exit(1); # programa i que aparegui un missatge d'error
}
my $tmp = ; # llegim la primera linia del format FASTA (l'identificador) i la guardem
@d = ; # llegim les linies de la sequencia, cada linia es una posicio de l'array
chomp $tmp; # treiem el canvi de linia que hi ha al final de l'identificador
$tmp=~/([\w>]+)/; # encaixem l'identificador amb l'expressio regular i
$id = $1; # el guardem a una variable singular (per evitar simbols estranys)
##### GUARDEM LA SEQUENCIA EN UN STRING #####
while ($c < scalar(@d)) { # llegim les linies de la sequencia i les guardem totes
chomp($d[$c]); # concatenades en una variable singular
$d[$c]=~/([acgtACGT]+)/; # eliminem qualsevol caracter que no sigui una base
$s = $s . $1;
$c = $c + 1;
}
$s = "\U$s"; # passem totes les bases a majuscules
close(SEQ); # tanquem el filehandle
return (\$id, \$s); # retornem en un vector la referencia de l'identificador i el vector en un string
}
# FUNCIO 2(mat): obre un argument que conte una matriu de pesos, l'emmagatzema en un hash
# i retorna la referencia al hash i el nombre de posicions que conte.
#
# Parametres requerits: la variable que emmagatzema el nom del fitxer de matriu de pesos que hem d'obrir.
sub mat {
my $ent= $_[0]; # variable que emmagatzema l'arxiu que hem d'obrir
my %mat; # hash on enregistrarem la matriu de pesos
my @n; # vector on enregistrarem els nucleotids
my $llegint=0; # variable indicadora de quan llegim la matriu
my $posicio=0; # posicio de la matriu que s'esta llegint
my $i=0; # comptador de la iteracio.
if (!open(MAT,"< $ent")) { # obrim la matriu de pesos emmagatzemada a l'argument
print " NO S'HA POGUT OBRIR EL FITXER DE LA MATRIU!!!\n"; # i, en cas que no es pugui obrir, fem saltar el
exit(1); # programa i que aparegui un missatge d'error
}
##### INTRODUIM LA MATRIU DE PESOS EN UN HASH #####
while () { # mentre llegim el fitxer de la matriu, si $llegint val 0
if ($llegint == 0 && m/\AP0\s+(\w+)\s+(\w+)\s+(\w+)\s+(\w+)/) { # i la linia que llegim comença amb P0
$n[0] = $1; # i encaixa amb l'expressio regular li direm que
$n[1] = $2; # emmagatzemi els elements que es troben en les diferents
$n[2] = $3; # columnes ($1, $2, etc) en les diferents posicions
$n[3] = $4; # del vector @n
$llegint = 1; # donem el valor de 1 a $llegint (es un flag)
}
if ($llegint == 1 && # si $llegint val 1 i encaixa la linia (des del
m/\A\d+\s+([\-\d\.]+)\s+([\-\d\.]+)\s+([\-\d\.]+)\s+([\-\d\.]+)/) { # principi) amb l'expressio regular, emmagatzemem
$mat{$n[0]}[$posicio] = $1; # en un hash els elements continguts en les
$mat{$n[1]}[$posicio] = $2; # diferents columnes, cada columna de les quals
$mat{$n[2]}[$posicio] = $3; # te una clau determinada, emmagatzemada al vector
$mat{$n[3]}[$posicio] = $4; # @n
$posicio = $posicio + 1; # cada cop que llegim una linia sumem +1 a $posicio
}
if ($llegint == 1 && m/\AXX\s+(\d+)/) { # si el flag val 1 i encaixa amb l'expressio regular
$mat{P} = $1; # donada, guardem el valor que hi ha a la primera
$llegint = 0; # columna en la clau P del hash
} # posem el flag a 0
}
close(MAT);
return (\$posicio , \%mat); # retornem el nombre d'elements que te la matriu i la matriu emmagatzemada en un hash
}
# FUNCIO 3(ORF): donada una sequencia en format fasta emmagatzemada en un vector, examina aquesta
# sequencia per a una determinada pauta de lectura en blocs de longitud determinada i retorna un
# hash on cada key es un array on s'emmagatzemen les posicions d'inici i de final d'un ORF
# (determinem un ORF com la sequencia que esta definida entre dos codons de STOP i, en el seu
# defecte, entre l'inici i el primer STOP o entre l'ultim STOP i el final de la sequencia).
#
# Parametres requerits: la sequencia en un string, la pauta de lectura.
sub ORF {
my $seq = $_[0]; # variable que emmagatzema la sequencia tota seguida
my $i = $_[1]; # inicialitzem el comptador de la iteracio amb el valor de la pauta de lectura que farem servir (0,1,2)
my $pos = 3; # emmagatzema el valor de la longitud de la finestra que mirem (en aquest cas 3)
my $s; # emmagatzema la sequencia tota seguida
my @stop; # guarda la posicio de l'ultim nucleotid de codo de stop trobat
my $finestra; # guarda el tros de sequencia que volem analitzar en la iteracio
$stop[@stop] = $i; # inicialitzem el vector que ens guarda la posicio dels codons de stop
# amb el valor de la posicio per on començarem a llegir la sequencia
# (determina l'inici del primer ORF)
##### BUSQUEM CODONS D'STOP #####
while ($i <= length ($seq) - $pos){
$finestra = substr ($seq, $i, $pos); # definim el tros de sequencia que explorem cada cop
if ($finestra =~m(TAA|TAG|TGA)){ # si fem un match d'un dels diferents codons de stop guardem
$stop[@stop] = $i + 2; # la ultima posicio del codo en un vector
}
$i = $i + 3; # augmentem tres el comptador (cada codo te tres elements)
}
$stop[@stop] = length ($seq); # assignem a l'ultim element del vector l'ultima posicio de la sequencia
my %results; # definim un hash en que emmagatzemarem els diferents ORF
my $count=0; # definim un comptador que sera les keys del hash
my @c; # definim el vector que emmagatzema els parells inici-final de ORF
$i = 0; # posem el comptador de la iteracio a 0
##### DEFINIM ORFs #####
while ($i < (@stop - 1)){ # fem un hash on cada key tindra un parell inici i
if ($i == 0){ # final de ORF
$results{$count}[0] = $stop[$i]; # resseguim el vector amb els codons de STOP,
} # la posició d'inici
# i la ultima posicio del vector
# guardem les diferents parelles dins d'un array
else{
$results{$count}[0] = $stop[$i]+1;
}
$results{$count}[1] = $stop[$i+1]; # introduim com a ultim element del vector la ultima
# posicio de la sequencia
$i++; # incrementem el comptador del while
$count++; # incrementem el comptador de les keys
}
return \%results; # retornem el hash amb els ORFs trobats
}
# FUNCIO 4(pes): fem una funcio que donada una matriu de pesos puntui tots els
# possibles senyals i retorni la posicio i la puntuacio de cada senyal.
#
# Parametres requerits (en aquest ordre): la sequencia en un string ($sequencia),
# el hash amb els ORFs, la matriu de pesos en un hash, el nombre de posicions
# de la matriu de pesos ($pos1) i el llindar que s'ha d'aplicar.L'ultim parametre
# de la funció indica si estem buscant donors (2) o acceptors(3) i la posicio en
# que s'ha d'emmagatzemar el hash resultant dins de l'array que esta emmagatzemat
# en cadascuna de les keys del hash %ORFs.
sub pes{
my $seq = $_[0]; # emmagatzemem la sequencia tota junta en un string
my %orf = %{$_[1]}; # guarda el hash amb els ORF trobats en un frame determinat
my %matriu = %{$_[2]}; # guarda el hash de la matriu de pesos que s'aplica
my $elements = $_[3]; # guarda el nombre d'elements que te la matriu que s'examina
my $llindar = $_[4]; # guarda la puntuacio minima que volem que tingui un senyal
my $tipus = $_[5]; # emmagatzema un escalar que indica la posicio en que hem
# d'emmagatzemar el hash resultant dins l'array de %ORF.
my $finestra; # emmagatzema el tros de seq que hem de mirar a cada volta,
# corresponent a un ORF
my @keysorf = keys(%orf); # guardem les claus del hash %orf en un array.
my $i = 0; # inicialitzem el comptador de la iteracio a zero.
my $h = 0; # inicialitzem el comptador de les keys del hash a zero.
########## RESSEGUEIX ORFs ##########
while ($i < @keysorf){
my %punt; # guarda la posicio i la puntuació dels senyals trobats en un ORF
my $count = 0; # comptador de les keys de %punt.
my @examinar = @{$orf{$h}};
my $inici = $examinar[0]; # definim les posicions d'inici i de final de
my $final = $examinar[1]; # cada ORF que volem analitzar.
if ($i == 0){ # si estem mirant el primer orf
if ($tipus == 2){ # estem buscant donors
$inici = $examinar[0]; # l'inici es igual a l'inici de l'ORF.
$final = $examinar[1] + 9; # el final es el final de l'ORF mes el nombre de posicions que te la matriu de pesos.
}
if ($tipus == 3){ # estem buscant acceptors
$inici = $examinar[0]; # l'inici es igual a l'inici de l'ORF.
$final = $examinar[1]; # el final es igual al final de l'ORF.
}
}
if ($i == (@keysorf - 1)){ # si estem mirant l'ultim ORF
if ($tipus == 2){ # estem buscant donors
$inici = $examinar[0]; # l'inici de la sequencia es igual al'inici de l'ORF.
$final = $examinar[1]; # el final de la sequencia es igual al final de l'ORF.
}
if ($tipus == 3){ # estem buscant acceptors
$inici = $examinar[0] - 21; # l'inici es l'inici de l'ORF menys el nombre de posicions que te la matriu de pesos.
$final = $examinar[1]; # el final de la sequencia es igual al final de l'ORF.
}
}
if ($i != 0 && $i != (@keysorf-1)) { # si mirem qualsevol altre ORF
if ($tipus == 2){ # estem buscant donors
$inici = $examinar[0]; # l'inici de la sequencia es igual al'inici de l'ORF.
$final = $examinar[1] + 9; # el final es el final de l'ORF mes el nombre de posicions que te la matriu de pesos.
}
if ($tipus == 3){ # estem buscant acceptors
$inici = $examinar[0] - 21; # l'inici es l'inici de l'ORF menys el nombre de posicions que te la matriu de pesos.
$final = $examinar[1]; # el final de la sequencia es igual al final de l'ORF.
}
}
my $llarg = $final - $inici +1; # definim la llargada de sequencia que utilitzarem per la cerca dels senyals
# (variara en funcio de la senyal i l'ORF).
$finestra = substr ($seq, $inici, $llarg); # definim el tros de sequencia que explorem cada cop i
my @vseq = split (//,$finestra); # el guardem en un vector
########## BUSCA SENYALS DINS D'UN ORF ##########
my $c = 0; # comptador del segon while a zero.
while ( $c < @vseq - $elements){ # analitzem la sequencia de cada ORF.
my $p = 0; # emmagatzema la puntuacio de la sequencia.
my $k = 0; # inicialitzem el comptador de posicio.
###### SUMA LA PUNTUACIO DE LES POSICIONS ######
while ( $k < $elements){ # analitzem un fragment de la sequencia de l'ORF de longitud.
$p = $p + $matriu{$vseq[$c + $k]}[$k]; # igual al nombre d'elements de la matriu de pesos i sumem.
$k = $k + 1; # les posicions dels diferents elements.
}
if ($p > $llindar){ # si la puntuació es superior al llindar donat
$punt{$count}[0] = $c + $matriu{P} + $inici; # emmagatzemem la posicio i la puntuacio de la
$punt{$count}[1] = $p; # senyal en un array
$count++; # incrementem el comptador de les keys del segon hash
}
$c++; # incrementem el comptador del segon while
}
$orf{$h}[$tipus]= \%punt; # guarda a cada volta la referencia al hash amb els resultats trobats
# en una posició determinada en funcio de les senyals que estem
# buscant dins el hash %orf
$i++; # incrementem el comptador del primer while
$h++; # incrementem el comptador de les keys del hash
}
return \%orf; # retornem el hash amb les senyals trobades.
}
# FUNCIO 5 (def): donades una serie de senyals donadores d'splicing i acceptores d'splicing
# fa els emparellaments de les diferents senyals acceptor-donor que es troben dins d'un mateix
# ORF (un acceptor amb tants donors com hi hagi al frame fins a un màxim determinat) i en suma
# les puntuacions.
# També defineix per cada ORF els emparellaments ATG-acceptor (determinen un exó inicial)
# i donor-STOP (defineixen un exo terminal)
# Retorna un array d'arrays on a cada subarray es un exo, del que retornem la posicio d'inici,
# la posicio de final i la seva puntuacio.
# Introduim a la funcio el hash ORF d'un frame determinat que conte un hash amb els diferents
# ORFs trobats al frame (que conte un array amb la posicio d'inici i final de l'ORF, un hash amb
# les senyals dels donors i la seva puntuació i un hash amb les senyals dels acceptors i la seva
# puntuació i el nombre d'emparellaments maxim que volem que faci.
#
# Parametres requerits: hash amb les senyals predites per cada ORF, nombre de prediccions,
# la sequencia, el frame, la longitud minima que volem per la prediccio.
sub def{
my %orf = %{$_[0]}; # guarda el hash amb les diferents senyals que hem de mirar.
my $max = $_[1]; # guarda el nombre maxim de prediccions que volem fer per cada senyal.
my $seq = $_[2]; # guarda la sequencia global que analitzem en una string.
my $frame = $_[3]; # guarda el frame en que estan les diferents senyals.
my $tamany = $_[4]; # guarda la longitud minima que volem que tingui un exo predit.
my $i = 0; # posem el comptador del primer while a 0.
my $results; # variable de treball que emmagatzema els resultats de les diferents funcions.
my @first; # guarda a cada posicio la referencia al array amb la informacio d'un exo first.
my @internal; # guarda a cada posicio la referencia al array amb la informacio d'un exo intern.
my @terminal; # guarda a cada posicio la referencia al array amb la informacio d'un exo final.
##### RESSEGUEIX ORFs #####
while ($i < keys(%orf)){ # mentre hi hagi ORFs
if ($orf{$i}[2]) { # si hi ha almenys un donor
my %atg; # guarda la primera posicio dels diferents ATG trobats
my @examinar = @{$orf{$i}};
my $inici = $examinar[0]; # definim les posicions d'inici i de final de
my $final = $examinar[1]; # cada ORF que volem analitzar
my $llarg = $final - $inici + 1; # definim la llargada de sequencia inclosa a cada ORF
my $finestra = substr ($seq, $inici, $llarg); # definim la sequencia que hi ha a cada ORF
##### BUSQUEM ATGs A CADA ORF #####
my $c = 0; # posem el comptador del segon while a 0;
my $pos = 0; # posem a cero el comptador de les keys de %atg.
while ($c <= length ($finestra) - 3){ # mirem la sequencia que hi ha dins l'ORF
my $busca = substr ($finestra, $c, 3); # per codons
if ($busca =~ m/ATG/g) { # si trobem un ATG emmagatzemem la posicio de
$atg{$pos}[0] = $c + $inici +1; # la primera base del codo la primera key del hash %atg.
$pos++; # incrementem el comptador de les keys del hash %atg.
}
$c = $c + 3; # augmentem tres el comptador (cada codo te tres elements)
}
##### EMPARELLEM ATG-DONOR #####
if (%atg && %{$orf{$i}[2]}){ # si s'ha trobat algun ATG en un frame que tambe te donors
##### CRIDEM LA FUNCIO PAR #####
$resultats = par(\%atg, \%{$orf{$i}[2]}, $max, $frame, $tamany); # assignem la referencia a l'array d'arrays
# que retorna la funcio par a $resultats.
my $n = 0; # inicialitzem a cero el comptador de la iteracio
my @mirem = @{$resultats}; # definim a cada volta l'array amb la informacio d'un exo que
# volem mirar.
while ($n < @mirem){ # resseguim els elements (referencies d'arrays) de @resultats
$first[@first] = $mirem[$n]; # guardem la referencia als diferents exons trobats en els
# diferents ORFs tots en el mateix array.
$n++; # incrementem el comptador de la primera iteracio
}
}
##### EMAPARELLEM ACEPTOR-DONOR #####
if ($orf{$i}[2] && $orf{$i}[3]) { # si hi ha almenys un donor i un acceptor per un orf determinat
##### CRIDEM LA FUNCIO PAR #####
$resultats = par (\%{$orf{$i}[3]}, \%{$orf{$i}[2]}, $max, $frame, $tamany); # assignem la referencia a l'array d'arrays
# que retorna la funcio par a $resultats.
my $n = 0; # inicialitzem el comptador de la iteracio a cero.
my @mirem = @{$resultats}; # definim a cada volta l'array amb la informacio
# d'un exo que volem mirar.
while ($n < @mirem){ # resseguim els elements (referencies d'arrays) de @mirem
$internal[@internal] = $mirem[$n]; # i assignem cadscuna de les subarrays en @internal (guardem
# les prediccions dels diferents ORFs en el mateix array).
$n++; # incrementem el comptador de la primera iteracio
}
}
##### EMPARELLEM DONOR-STOP #####
if (%{$orf{$i}[3]}){ # mentre hi hagi algun acceptor.
my %donors = %{$orf{$i}[3]}; # assignem a cada volta el hash que hi ha als donors.
my @n = keys (%{$orf{$i}[3]}); # guardem les keys del hash en un vector.
my $k = scalar (@n); # guardem el nombre d'elements que te el hash.
my $final = $orf {$i}[1]; # guarda la posició de final d'exo (que es la de final d'Orf).
my $c = 0; # posem el comptador de la iteracio a cero.
my @resultats; # guarda la referencia als diferents exons trobats
while ($c < $max){ # mentre no hagim fet el nombre d'emparellaments maxim
if ($donors{$k}){ # si hi ha algun donor
my @exo; # definim el vector exo
$exo[0] = $donors{$k}[0]; # guardem la posicio d'inici de l'exo a la primera posicio del vector
$exo[1] = $final + 1; # guardem la posicio de final d'exo a la segona posicio del vector
$exo[2] = $donors{$k}[1]; # guardem la puntuacio de l'exo a la tercera posicio del vector
$exo[3] = $frame; # guardem el frame en que s'ha predit l'exo.
$exo[4] = "Terminal"; # guardem l'identificador de tipus d'exo.
if (($exo[1] - $exo[0] + 1) > $tamany){ # si la longitud de la prediccio es superior a la minima definida
$terminal[@terminal] = \@exo; # guardem la referencia a cadascun dels exons terminals predits en
} # el tots els ORFs en el vector @terminal.
}
$c++; # incrementem el comptador del hash
$k--; # disminuim el comptador de les claus dels ORFs
}
}
}
$i++; # incrementem el comptador de les keys del hash %orf
}
return (\@first, \@internal, \@terminal); # retornem les referencies als arrays pels diferents tipus d'exons predits
}
# FUNCIO 6 (par): realitza els diferents emparellaments que es puguin fer entre dos tipus de
# senyals definides cadascuna d'elles en un hash i retorna l'emparellament i la puntuació de
# l'exo trobat en un hash, on cadascuna de les claus conte un array amb l'inici, el final i
# la puntuació del exo predit.
#
# Parametres requerits: hash amb les senyals d'inici d'exo, hash amb les senyals de final, nombre
# de prediccions maximes per senyal, frame en que hem trobat l'exo i longitud minima de l'exo.
sub par{
my %seni = %{$_[0]}; # guarda el hash amb les senyals que determinen l'inici d'exo.
my %senf = %{$_[1]}; # guarda el hash amb les senyals que determinen final d'exo.
my $max = $_[2]; # guarda el nombre màxim d'emparellaments que volem fer per cada senyal d'inici.
my $frame = $_[3]; # guarda el frame en que s'ha predit les senyals
my $tamany = $_[4]; # longitud minima que volem que tingui l'exo predit
my $f = 0; # posem el comptador del primer while a zero.
my @resultats; # guarda els arrays amb la informacio de un exo en un array.
##### MIREM LES SENYALS D'INICI #####
while ($f < keys (%seni)){ # mentre hi hagi senyals d'inici
my @inici = @{$seni{$f}}; # definim l'array que conte la informacio de la senyal d'inici
my $k = 0; # posem el while del segon comptador a zero.
my $count = 0; # posem el comptador dels emparellaments a zero.
##### MIREM LES SENYALS DE FINAL #####
while ($k < keys (%senf)){ # examinen el hash %senf mentre hi hagi senyals
my @final = @{$senf{$k}}; # definim a cada volta l'array que hi ha dins de
# la key del hash %senf
if ($seni{$f}[0] < $final[0]){ # si la senyal d'inici d'exo esta abans que la de final d'exo
my @exo; # definim un array que guardara la informacio del exo
$exo[0] = $inici[0]; # guardem a la primera posicio la senyal d'inici de l'exo
$exo[1] = $final[0]; # guardem a la segona posicio la senyal de final de l'exo
if ($inici[1]){ # si existeix la posicio $inici[1](l'inici ve definit per un acceptor)
$exo[2] = $inici[1] + $final[1]; # la tercera posicio es la suma de les puntuacions de la senyal
$exo[4] = "Internal"; # d'inici i la de final (donor).
}
else{ # si no existeix la posicio $inici[1](l'inici ve definit per un ATG
$exo[2] = $final[1]; # la tercera posicio del vector es igual a la puntuacio de la senyal
$exo[4] = "First "; # de final (donor).
}
$exo[3] = $frame;
if ($count < $max && ($exo[1] - $exo[0] + 1) > $tamany){ # guardem els resultats dels exons predit fins a un nombre
$resultats[@resultats] = \@exo; # maxim determinat per $max en @resultats.
$count++; # incrementem el comptador d'@resultats
}
}
$k++; # incrementem el comptador del while.
}
$f++; # incrementem el comptador de les keys de les senyals de final d'exo
}
return \@resultats; # retornem un array amb les referencies als arrays dels diferents exons trobats.
}
# FUNCIO 7 (bcod): donat un fitxer que conte una matriu de biaix codificant, l'obre i
# i l'emmagatzema en un hash. Retorna el hash que conte la matriu de codons.
sub bcod {
my $fitxercod = $_[0]; # variable que emmagatzema les distribucions de probabilitat dels codons
my %pcodons; # guarda el hash amb el biaix codificant
if (!open(PROPCODONS,"< $fitxercod")) { # comprova que s'obri el fitxer que conte la matriu de codons i sino
print "impossible obrir $fitxercod\n"; # fa saltar el programa i un missatge d'error.
exit(1);
}
##### LLEGIM EL FITXER I L'EMMAGATZEMEM EN UN HASH #####
while () { # obrim el filehandle
$_ =~m/(\w+) (\d+\.\d+)/; # expressio regular que encaixa amb dues
$pcodons{$1} = $2; # paraules la primera de les quals es alfabetica $1 i
# la segona un numero amb decimals $2
}
close(PROPCODONS); # tanquem el filehandle
return \%pcodons; # retorna el hash amb la matriu introduida
}
# FUNCIO 8 (biaix): calcula la rao de versemblança dels possibles exons en una pauta de
# lectura donada i retorna el valor d'aquesta rao en escala logaritmica.
#
# Parametres requerits: sequencia en una string, el vector amb la informacio de l'exo
# (inici, final i puntuacio) i el frame. Retorna el mateix vector amb la puntuacio
# del biaix codificant inclosa nomes en cas de que sigui superior al llindar
sub biaix {
my $seq = $_[0]; # emmagatzemem la sequencia emmagatzmada en un string.
my %codons = %{$_[1]}; # emmagatzema el hash amb el biaix codificant.
my @prediccions = @{$_[2]}; # guarda l'array que conte els arrays amb les prediccions dels diferents tipus de
# gens trobats per un ORF determinat.
my $llindar = $_[3]; # guarda la puntucio minima global que volem que tingui una prediccio.
my $i = 0; # inicialitzem el comptador del while a zero.
my @resultats; # guarda l'array que conte els arrays dels diferents exons que puntuin per sobre del llindar
##### RESSEGUIM ELS ELEMENTS DE L'ARRAY @PREDICCIONS #####
while ($i < 3){ # mentre hi hagi un tipus d'exo (nomes mirem inicials,interns i terminals)
my @mirem = @{$prediccions[$i]}; # definim a cada volta l'array que te les prediccions per
# un tipus d'exo determinat.
my $n = 0; # inicialitzem el comptador de la segona iteracio a zero.
##### RESSEGUIM ELS DIFERENTS EXONS QUE HI HA DE CADA TIPUS #####
while ($n < @mirem){ # mentre hi hagi prediccions d'un tipus d'exo determinat
my @exo = @{$mirem[$n]}; # definim l'array amb cadascun dels exons
my $longexo = ($exo[1]-1) - ($exo[0] - 1) + 1; # trobem la longitud de l'exo
my $seqexo = substr($seq,$exo[0] - 1,$longexo); # prenem el tros de sequencia corresponent a l'exo
my $frame = ($exo[0] - 1) - $exo[3]; # determinem la longitud de la sequencia
##### DETERMINEM EL FRAME DE L'EXO #####
if ($frame % 3 == 0){ # si la posicio d'inici del nostre exo es el de la primera
$exo[3] = 0; # base del codo, el frame de l'exo sera 0.
}
if ($frame % 3 == 1){ # si la posicio d'inici del nostre exo es el de la segona
$exo[3] = 2; # base del codo, el frame de l'exo sera 2.
}
if ($frame % 3 == 2){ # si la posicio d'inici del nostre exo es el de la primera
$exo[3] = 1; # base del codo, el frame de l'exo sera 1.
}
if (length ($seqexo) > 2){ # si la longitud es suficient com per aplicar el hash de codons
my $seq_exo = substr($seqexo, $exo[3]); # definim la longitud de sequencia a examinar (en funcio del frame)
my @codo = ($seq_exo =~ m/.../g); # guardem la sequencia en un array dividida per codon.
my $r = 0; # inicialitzem a zero la variable que emmagatzema la puntuacio de cada prediccio.
my $c = 0; # inicialitzem a zero el comptador de la segona iteracio.
##### RESSEGUIM LA SEQUENCIA DE L'EXO PREDIT #####
while ($c < scalar(@codo)) { # mentre hi hagi codons
$r = $r + log($codons{$codo[$c]}) - log(1/64); # la puntuacio es la suma de la que trobem aplicant
# el biaix codificant sobre cada codo de la sequencia
$c++; # incrementem el comptador del tercer while.
}
$exo[2]= $exo[2]+$r; # la puntuacio de l'exo es la puntuacio del biaix
# codificant mes la de les senyals d'splicing.
if ($exo[2] > $llindar){ # si la puntuacio es superior a la que hem determinat
# com a minima per acceptar una prediccio.
$prediccions[$i][$n][2] = $exo[2]; # emmagatzemem la nova puntuacio a la posicio que contenia
# la puntuacio de les senyals a @prediccions.
# guardem el contingut del vector que hi ha a
$resultats[$i][$n][0] = $exo[0]; # @prediccions i guardem els diferents parametres que conte
$resultats[$i][$n][1] = $exo[1]; # al vector @resultats (inici, final, puntuacio i frame)
$resultats[$i][$n][2] = $exo[2];
$resultats[$i][$n][3] = $exo[3];
$resultats[$i][$n][4] = $exo[4];
}
}
$n++; # incrementem el comptador del segon while.
}
$i++; # incrementem el comptador del primer while
}
return \@resultats; # retornem la referencia a l'array amb les prediccions fetes
# amb la nova puntuacio
}
# FUNCIO 9 (junta): donat un array d'array d'arrays, emmagatzema tots els elements
# en un array d'arrays (treiem un nivell de complexitat).
#
# Parametres requerits: l'array que volguem arreglar.
sub junta {
my @array = @{$_[0]}; # emmagatzema el vector que volem arreglar
my $i = 0; # inicialitzem el comptador de la iteracio a zero.
my @resultats; # emmagatzema l'array amb un nivell menys de complexitat.
while ($i < @array){ # mentre hi hagi elements del primer nivell
if ($array[$i]) { # si existeix un element concret
my $c = 0; # inicialitzem el comptador de la iteracio a zero.
my @primer = @{$array[$i]}; # definim a cada volta l'array amb les referencies
while ($c < @primer){ # mentre hi hagi elements a cada subarray
if ($primer[$c]){ # si existeix un element exo en una posicio concreta
my @segon = @{$primer[$c]};
$resultats[@resultats] = \@segon; # guarda tots els elements seguits en @resultats
}
$c++; # incrementem el comptador del segon while
}
}
$i++; # incrementem el comptador del primer while
}
return \@resultats; # retornem l'array arreglat.
}
###################### FI DEL PROGRAMA EXON-finding #####################################
