#!/usr/bin/perl
use strict;


### Programa escrit amb llenguatge PERL 

#################### CERCAEXONS  #####################################
########## Treball de bioinformatica. Universitat Pompeu Fabra #######
######### Alexandre  Gomez Garrido i Alba Vidal Torres ###############
################ Barcelona, marc 2005 ################################
######################################################################


################# PARAMETRES DE ENTRADA ##############################
### El programa necessita la introduccio dels seguents parametres:
### Nom del fitxer fasta
### Nom del fitxer de la matriu d'acceptors 
### Nom del fitxer de la matriu de donadors 
### Nom del fitxer de la taula d'us de codons 
### Valor del llindar a partir del qual considerem un lloc de splicing candidat com a lloc de splicing real
### Valor del llindar a partir del qual considerem un exo predit com a possible exo real
### Exemple d'execucio del programa: ./cercaexons.pl seq.fa matriu_acceptors.txt matriu_donadors.txt taula_biaix.txt 1 0 
#####################################################################




### Comprovem que s'hagin introduit els 7 parametres d'entrada 

if (scalar(@ARGV) != 6) {			
	print "\nError:  Has d'introduir els parametres en aquest ordre: \n\tFitxer de la sequencia en format FASTA.\n\tMatriu de pesos de llocs de splicing acceptors.\n\tMatriu de pesos de llocs de splicing donadors.\n\tTaula d'us de codons.\n\tLLindar pels llocs de splicing.\n\tLLindar pels exons.\n\n";
	exit(-1);
}


### DECLARACIO DE LES VARIABLES 

### Enregistrem els parametres d'entrada en variables mes entenedores 

my $fasta = $ARGV[0];			### Conte el nom del fitxer fasta 
my $acceptors = $ARGV[1];		### Conte el nom del fitxer de la matriu d'acceptors 
my $donadors = $ARGV[2];		### Conte el nom del fitxer de la matriu de donadors 
my $codons = $ARGV[3];			### Conte el nom del fitxer de la taula d'us de codons 
my $llindar = $ARGV[4];			### Conte el valor del llindar a partir del qual considerem un lloc de splicing candidat com a lloc de splicing real	
my $llindar_exons = $ARGV[5];		### Conte el valor del llindar a partir del qual considerem un exo predit com a possible exo real 

my $id;					### Enregistrarem l'identificador de la sequencia entrada 
my @s;					### Enregistrarem la sequencia 
my @s_reverse;				### Enregistrarem la sequencia amb strand negatiu 
my %matriu_acceptors;			### Hash on enregistrarem la matriu de pesos dels llocs de splicing acceptors 
my %matriu_donadors;			### Hash on enregistrarem la matriu de pesos dels llocs de splicing donadors 
my %taula_codons;			### Hash on enregistrarem la taula d'us de codons 
my @resultat_pos;			### Matriu on enregistrarem els exons predits per l'strand positiu
my @resultat_neg;			### Matriu on enregistrarem els exons predits per l'strand negatiu
my @resultat_pos_bons;			### Matriu on enregistrarem els exons predits per l'strand positiu que hagin superat el llindar
my @resultat_neg_bons;			### Matriu on enregistrarem els exons predits per l'strand negatiu que hagin superat el llindar
my @resultats_units;			### Matriu que contindra els exons predits


### LECTURA DELS FITXERS D'ENTRADA I ENREGISTRACIO 

($id, @s) = LLEGIR_FASTA($fasta);		### Al  vector @s li assignem la sequencia 
%matriu_acceptors = LLEGIR_MATRIU($acceptors);
%matriu_donadors = LLEGIR_MATRIU($donadors);
%taula_codons = LLEGIR_USCODONS($codons);


@s_reverse = reverse(@s);			### Enregistrem la sequencia amb strand negatiu en un nou vector 


### BUSQUEDA D'EXONS PER ELS DOS STRANDS

@resultat_pos = EXONS(\@s, \%matriu_acceptors, \%matriu_donadors, \%taula_codons, $llindar, "+");

@resultat_neg = EXONS(\@s_reverse, \%matriu_acceptors, \%matriu_donadors, \%taula_codons, $llindar, "-");

### TRADUCCIO DE LES POSICIONS DE L'STRAND NEGATIU, PER LES POSICIONS REALS EN LA SEQUENCIA ORIGINAL (STRAND POSITIU)

@resultat_neg = POSICIONS_REALS_STRANDNEG (\@resultat_neg, scalar(@s)); 

### SELECCIO DELS RESULTATS QUE SUPEREN EL LLINDAR ESTABLERT PELS EXONS

@resultat_pos_bons = EXONS_BONS(\@resultat_pos, $llindar_exons); 
@resultat_neg_bons = EXONS_BONS(\@resultat_neg, $llindar_exons);


### UNIO DE LES DUES MATRIUS DE RESULTATS EN UNA MATRIU FINAL

@resultats_units = UNIR_RESULTATS(\@resultat_pos_bons, \@resultat_neg_bons);


### ORDENACIO DE LA MATRIU DE MENOR A MAJOR POSCIO DE NUCLEOTID I EN CAS D'EMPAT, DE MAJOR A MENOR PUNTUACIO

@resultats_units = sort{${$a}[0] <=> ${$b}[0] ? ${$a}[0] <=> ${$b}[0] :
${$a}[1] <=> ${$b}[1] ?${$a}[1] <=> ${$b}[1] :
 ${$b}[2] <=> ${$a}[2]}(@resultats_units);

### ESCRIPUTURA DEL FITXER DE SORTIDA

FORMAT_GFF($id,\@resultats_units);
 

 
########################## EXONS #####################################
### FUNCIO PRINCIPAL DEL PROGRAMA. L'EXECUTEM PER LES DUES STRANDS ###
###################################################################### 
### Aquesta funcio permet trobar tots els possibles exons donada una sequencia
### Parametres d'entrada:
######################### Adreca de la sequencia
######################### Adreca del hash de la matriu de acceptors
######################### Adreca del hash de la matriu de donadors
######################### Adreca del hash de la taula d'us de codons
######################### Valor llindar a partir del qual considerem un lloc de splicing candidat com a lloc de splicing real
######################### Strand
### Retorna: Matriu de resultats que conte els exons predits
#####################################################################


sub EXONS {
	
	
### DECLARACIO VARIABLES 
	
	my @s = @{$_[0]};
	my %matriu_acceptors = %{$_[1]};
	my %matriu_donadors = %{$_[2]};
	my %taula_codons  = %{$_[3]};
	my $llindar = $_[4];
	my $strand = $_[5];
	
### Vectors de sequencia en la pauta corresponent, on a cada posicio del vector enregistrem un codo
	my @s0;				
	my @s1;		
	my @s2;
				
### Matrius on a cada fila contenen una sequencia de nucleotids on els tres darrers nucleotids son un codo stop (en la pauta corresponent)
	my @matriu0;	
	my @matriu1;
	my @matriu2;
	
### Matrius on enregistrarem els resultats per tots els exons predits
	my @resultat0;
	my @resultat1;
	my @resultat2;
	
### Matriu que conte tots els resultats. Es el vector que retornara aquesta funcio
	my @resultatJ;
	
### Generem 3 vectors on a cada posicio contenen un codo, segons la pauta de lectura corresponent 
	@s0 = PAUTES_LECTURA(\@s, 0);
	@s1 = PAUTES_LECTURA(\@s, 1);
	@s2 = PAUTES_LECTURA(\@s, 2);
	
### Generem 3 matrius on enregistrem la sequencia separada per codons STOP
	@matriu0 = MATRIU_STOPCODONS(\@s0,\@s,0);
	@matriu1 = MATRIU_STOPCODONS(\@s1,\@s,1);
	@matriu2 = MATRIU_STOPCODONS(\@s2,\@s,2);
		
### Predim exons i els enregistrem a les matrius de resultats
	@resultat0 = RESULTAT(\@matriu0, \%matriu_acceptors, \%matriu_donadors, $matriu_acceptors{P}, $matriu_donadors{P}, \%taula_codons, 0, $llindar, $strand, \@s);
	@resultat1 = RESULTAT(\@matriu1,  \%matriu_acceptors, \%matriu_donadors, $matriu_acceptors{P}, $matriu_donadors{P}, \%taula_codons, 1, $llindar, $strand, \@s);
	@resultat2 = RESULTAT(\@matriu2,  \%matriu_acceptors, \%matriu_donadors, $matriu_acceptors{P}, $matriu_donadors{P}, \%taula_codons, 2, $llindar, $strand, \@s);
	
### Matriu que conte tots els exons predits
	@resultatJ = RESULTATS_JUNTS(\@resultat0, \@resultat1, \@resultat2);
	
	return(@resultatJ);
}
     
     
#####################################################################
####################### FUNCIONS ####################################
#####################################################################

######################## LLEGIR_FASTA ############################### 
### Aquesta funcio obre i llegeix un fitxer en format FASTA i enregistra
### l'identificador de la sequencia i en un vector on a cada posicio hi 
### ha un nucleotid de la sequencia.
### Parametres d'entrada:
######################### Nom del fitxer FASTA
### Retorna: Vector de la sequencia i l'identificador
#####################################################################


sub LLEGIR_FASTA {

    if (!open(FASTA,"<$_[0]")) {	### Obrim el fitxer FASTA i sino es pot obrir donem un error
		print "Error: No es pot obrir el fitxer FASTA $_[0]";
		exit(-1);
	}		
    my $nom= <FASTA>;			### A $nom li assignem la primera linia del fitxer que conte l'identificador de la sequencia
        chomp($nom);
    if ($nom =~ m/\>([\d\w]+)/) {	### A l'identificador no agafem el simbol >
    	$nom = $1;
	}
    else {				### Si la primera linia del fitxer no compleix l'expressio regular, considerem que no esta en formar FASTA i imprimim en pantalla un error
    	print "error: el fitxer no es en format FASTA\n";
	exit(-1);
	} 
    
    my $seq; 				### Enregistrem la sequencia tota seguida en un variable singular
    my @s; 				### Enregistrarem la sequecia en un vector
    
  while(<FASTA>) {
	if ($_ =~ m/([ACTGactg]+)/) {
	    $seq = $seq . $1;                     ### Anem concatenant les linies de la sequencia
	}
    }
    
    $seq = uc($seq); 			### Passem tots els nucleotids a majuscules
    @s=split(//,$seq); 			### Separem els nucleotids un per un i els coloquem a una posicio del vector @s
	
    close(FASTA);			### Tanquem el fitxer FASTA
    return($nom, @s);			### Retornem l'identificador i el vector de la sequencia
}



####################### LLEGIR_MATRIU ############################### 
### Aquesta funcio obre i llegeix un fitxer que conte una matriu de pesos i l'enregistra en un hash. Cada clau correspon a un nucleotid. A cada nucleotid li assignem un vector que conte les frequencies en que trobem aquell nucleotid en cada posicio
### Parametres d'entrada:
######################### Nom del fitxer que conte la matriu de pesos
### Retorna: Hash amb les dades de la matriu de pesos
#####################################################################


sub LLEGIR_MATRIU {
    my %matriu_pesos; 			  ### Declarem el hash
    if (!open(MATRIU,"<$_[0]")){	  ### Obrim el fitxer i sino pot obrir-lo ens dona un error i s'acaba l'execucio
    	print "Error: No es pot obrir el fitxer de la matriu de pesos $_[0]";
	exit(-1);
    }

    my $temp=<MATRIU>;			  ### Llegim la primera linia i la coloquem a una variable singular
    my @nucleotids = split(/\s+/, $temp); ### Separem per espais per la primera linia, per tant el vector conte el nom dels nucleotids i el nom P0
    my $posicio;			  ### Enregistrem la posicio de la matriu en la que ens trobem
	
		
    while(<MATRIU>){                      
     	if ($_ =~ m/\A\d+\s+([\-\d\.]+)\s+([\-\d\.]+)\s+([\-\d\.]+)\s+([\-\d\.]+)/) { 

### Enregistrarem les frequencies de cada nucleotid en aquella posicio	
		
		$matriu_pesos{$nucleotids[1]}[$posicio] = $1;   
		$matriu_pesos{$nucleotids[2]}[$posicio] = $2;		
		$matriu_pesos{$nucleotids[3]}[$posicio] = $3;		
		$matriu_pesos{$nucleotids[4]}[$posicio] = $4;		
				
		$posicio = $posicio + 1;
        }
	$matriu_pesos{tamany} = $posicio;
	if ($_=~ m/\AXX\s+(\d+)/) {	   ### Complira aquesta expressio regular la ultima linia. Conte el nucleotid a partir del que comenca o acaba de codificar
		$matriu_pesos{P} =$1;	   ### Creem una nova clau que conte el numero del nucleotid
		}            
	}
          
    close(MATRIU);
    return(%matriu_pesos);
}

####################### LLEGIR_USCODONS ############################ 
### Aquesta funcio obre i llegeix un fitxer que conte la frequencia d'us de codons. Cada clau correspon un codo al que li assignem la seva frequencia
### Parametres d'entrada:
######################### Nom del fitxer que conte la taula d'us de codons
### Retorna: Hash amb les dades de la taula de codons
####################################################################

sub LLEGIR_USCODONS {

    my %us_codons;			### Declarem el hash
    my @linea;				### Declarem el vector on enregistrarem les dades de cada linia
    
    if (!open(CODONS,"<$_[0]")) {	### Obrim el fitxer que conte la taula d'us de codons i sino pot ens dona un error
    	print "Error: No es pot obrir el fitxer $_[0] amb la taula d'us de codons.";
	exit(-1);	
    }
    while (<CODONS>) {
        @linea = split(/\s+/,$_);  	### Separem cada linia per espais, a la primera posicio del vector tindrem el codo i a la segona la frequencia
        $us_codons{$linea[0]} = $linea[1]; ### La primera posicio del vector es la clau del hash, i la segona el valor que li assignem
        }

    close(CODONS);
    return(%us_codons);
}


######################### PAUTES_LECTURA ########################## 
### Aquesta funcio crea un vector on a cada posicio conte un codo, corresponent a una pauta de lectura determinada de la sequencia
### Parametres d'entrada:
######################### Adreca del vector de la sequencia, a cada posicio hi ha un nucleotid
######################### Pauta de lectura per la que volem crear la divisio de codons
### Retorna: Vector de codons segons la pauta introduida
####################################################################


sub PAUTES_LECTURA {

	my @s = @{$_[0]};		### Vector amb la sequencia
	my $i = $_[1];			### Pauta de lectura
	my @sp;				### Vector on enregistrarem els codons
	my $c = 0;
	
	while ($i < scalar(@s)-2) {			### Des del valor de la pauta de lectura fins a l'ultim triplet possible
		$sp[$c] = $s[$i].$s[$i+1].$s[$i+2];	### Concatenem 3 nucleotids per fer un codo
		$i = $i + 3;				### Aquest comptador va de 3 en 3, seguint l'esquema de la linia de adalt
		$c = $c + 1;				### Comptador de les posicions del vector de codons
		}

	return(@sp);
}


#################### MATRIU_STOPCODONS ############################ 
### Aquesta funcio crea una matriu. Cada fila es un vector de nucleotids en que els 3 ultims nucleotids formen un codo STOP, excepte l'ultima fila on tenim el tros de sequencia restant (pot o no acabar en un codo STOP). Aixi tenim enregistrada tota la sequencia originial pero separada per codons STOP.
### Parametres d'entrada:
######################### Adreca del vector que conte els codons corresponents a una pauta de lectura
######################### Adreca del vector de la sequencia, a cada posicio hi ha un nucleotid
######################### Pauta de lectura en la que hem fet la separacio de codons
### Retorna: Matriu de la sequencia separada per codons STOP.
####################################################################


sub MATRIU_STOPCODONS {

	my @scodons = @{$_[0]};			### Vector amb els codons
	my @s = @{$_[1]};			### Vector amb la sequencia nucleotid a nucleotid
	     
	my $pauta = $_[2];			### Pauta de lectura
	
	my $i = 0; 				### Comptador que recorrera tot el vector de codons
	my $inicial = 0;			### Nucleotid inicial d'una fila de la matriu
	my $final;				### Nucleotid inicial de la seguent fila de la matriu
	my @matriu;				### Matriu on enregistrem els resultats
	my $fila = 0;
	my $columna;
	my $n;					### Comptador que recorre el vector de nucleotids des de $inicial fins a $final
	
	while ($i<scalar(@scodons)) {
		if (($scodons[$i] =~ m/(TAA|TAG|TGA)/) || ($i == (scalar(@scodons)-1))) {  ### Si el codo es un codo stop o es l'ultim codo dela sequencia, enregistrem el fragment  
			$columna = 0;
			$n = $inicial;
			$final = (($i + 1)*3 + $pauta); ### $final pren el valor de la posicio nucleotid de despres del codo stop trobat
			while ($n < $final) {
				$matriu[$fila][$columna] = $s[$n]; ### Anem enregistrant el fragment de la sequencia en la matriu
				$n = $n + 1;
				$columna = $columna + 1;
				}
			if ($i == (scalar(@scodons)-1)) {	
				$matriu[$fila][$columna] = $s[$n];  #perque tambe ens escrigui l'ultim nucleotid
				} 
			$fila = $fila + 1;
			$inicial = $final;
		}
		$i = $i + 1;
	}
	return(@matriu);
}


#################### RESULTAT ################################ 
### Aquesta funcio crida varies funcions. El que finalment obtenim es una matriu amb les dades dels exons predits, a partir d'haver separat la sequencia per codons STOP. Aixo ho aconseguim a partir de recorrer la matriu (amb la sequencia separada per codons STOP) fila per fila, i buscar-hi tots els exons possibles, tenint en compte els llocs de splicing, els codons ATG, els codons STOP, i el biaix codificant dels possibles exons. 
### Parametres d'entrada:
######################### Adreca de la matriu que conte la sequencia separada per codons STOP
######################### Adreca del hash on hi ha les dades de la matriu de pesos de llocs de splicing acceptors
######################### Adreca del hash on hi ha les dades de la matriu de pesos de llocs de splicing donadors
######################### Posicio del lloc de splicing en la que comenca a codificar per un exo
######################### Posicio del lloc de splicing en la que acaba de codificar per un exo
######################### Adreca del hash on hi ha les dades de la taula d'us de codons
######################### Pauta de lectura per la que estem executant la funcio
######################### Valor del llindar a partir del qual considerem un lloc de splicing candidat com a lloc de splicing real
######################### Strand pel que estem executant les funcions
### Retorna: Matriu amb les dades dels exons predits
####################################################################

sub RESULTAT {

	my @s = @{$_[0]};
	my %hash_acceptors = %{$_[1]};
	my %hash_donadors = %{$_[2]};
	my $pa = $_[3];
	my $pd = $_[4];
	my %hash_biaix = %{$_[5]};
	my $pauta = $_[6];
	my $llindar = $_[7];
	my $strand = $_[8];
	my @seq = @{$_[9]};
	
	my @splicing_acceptors;		### Vector on enregistrarem les puntuacions dels candidats a llocs de splicing acceptors
	my @splicing_donadors;		### Vector on enregistrarem les puntuacions dels candidats a llocs de splicing donadors
	my @taula_acceptors;		### Matriu on enregistrarem les dades dels possibles llocs de splicing acceptors
	my @taula_donadors;		### Matriu on enregistrarem les dades dels possibles llocs de splicing donadors
	my @taula_ATG;			### Matriu on enregistrarem les posicions dels codons d'inici
	my @taula_STOP;			### Matriu on enregistrarem les posicions dels codons STOP
	my @combinacio_interns;		### Matriu on enregistrem les combinacions que donen lloc a possibles exons interns
	my @combinacio_inici;		### Matriu on enregistrem les combinacions que donen lloc a possibles exons inicials
	my @combinacio_term;		### Matriu on enregistrem les combinacions que donen lloc a possibles exons terminals
	my @combinacio_single;		### Matriu on enregistrem les combinacions que donen lloc a possibles exons simples (estan formats per un unic exo) 
	my @combinacio_3;		### Matriu on enregistrem tres de les combinacions anteriors (aixo ho fem per poder reutilitzar les funcions que ja tenim fetes per unir matrius)
	my @combinacio;			### Matriu on enregistrem totes les combinacions anteriors
	my @resultat_fila;		### Matriu on tenim les dades dels exons predits per una fila
	my @resultat;			### Matriu amb les dades dels exons predits
	
	my $i = 0;
	my $c = 0;
	my $j;
	my $cont;
	my $p;
	
	while ($i < scalar(@s)) {	### Recorrem totes les files de la matriu de la sequencia separada per codons STOP
		$cont = 0;
		$p = 0;
		while ($cont < $i) {			### Mentres no arribem a la darrera fila
		$p = $p + scalar(@{$s[$cont]});		### Anirem sumant tots els nucleotids d'aquella fila
		$cont = $cont + 1;
		}
		
		@splicing_acceptors = FILA_SPLICING(\@seq, \%hash_acceptors, $hash_acceptors{tamany}, $pa, $p, scalar(@{$s[$i]}));
		@taula_acceptors = SOBRE_LLINDAR(\@splicing_acceptors, $llindar); ### Escollim els llocs de splicing acceptors que superen el llindar
		
		@splicing_donadors = FILA_SPLICING(\@seq, \%hash_donadors, $hash_donadors{tamany}, $pd, $p, scalar(@{$s[$i]}));
		@taula_donadors = SOBRE_LLINDAR(\@splicing_donadors, $llindar);   ### Escollim els llocs de splicing acceptors que superen el llindar

		@taula_ATG = MATRIU_ATG($s[$i], $pauta, $i);	### Busquem tots els codons ATG
		@taula_STOP = MATRIU_STOP($s[$i]); 		### Enregistrem el codo STOP (nomes en podem tenir un per fila)

### Fem totes les combinacions requerides, per poder predir tots els possibles exons.
		@combinacio_interns = COMBINACIONS(\@taula_acceptors, \@taula_donadors, "internal"); 
		@combinacio_inici = COMBINACIONS(\@taula_ATG, \@taula_donadors, "first");
		@combinacio_term = COMBINACIONS(\@taula_acceptors,\@taula_STOP, "terminal");
		@combinacio_single = COMBINACIONS(\@taula_ATG,\@taula_STOP, "single");
		@combinacio_3 = RESULTATS_JUNTS(\@combinacio_inici, \@combinacio_interns, \@combinacio_term);
		@combinacio = UNIR_RESULTATS(\@combinacio_3, \@combinacio_single);
	
		@resultat_fila = BIAIX_CODIFICANT(\@combinacio, \%hash_biaix, $s[$i], $pauta, $i); ### Calculem el valor de biaix codificant per totes les sequencies seleccionades com a possibles exons
		@resultat_fila = POSICIONS_REALS(\@resultat_fila,\@s, $i, $strand); ### Traduim les posicions dels nucleotids a les posicions reals, ja que les tenim enregistrades respecte la seva fila

### Ara ajuntem els resultats que anem obtenint per les diferents files de la matriu que conte la sequencia separada per STOP codons.
		$j = 0;
		while ($j < scalar(@resultat_fila)) {
			$resultat[$c] = $resultat_fila[$j];
			$j = $j + 1;
			$c = $c + 1;
		}
		$i = $i + 1;
	}

	return(@resultat);
}

#################### FILA_SPLICING ################################ 
### Aquesta funcio crea un vector que conte les puntuacions de la fila. Cada posicio del vector conte la puntuacio del logaritme de likelihood de que en aquell nucleotid comenci un lloc de splicing. D'aquesta manera recorrem tots els nucleotids i calculem tots els possibles llocs de splicing.
### Parametres d'entrada:
######################### Adreca de la fila de la matriu que conte la sequencia que volem puntuar
######################### Adreca del hash on hi ha les dades de la matriu de pesos
######################### LLargada del lloc de splicing
### Retorna: Matriu de puntuacions
####################################################################
	
	
sub FILA_SPLICING {
	my @seq = @{$_[0]};			### Vector que conte una fila
	my %matriu = %{$_[1]};			### Hash de la matriu de pesos
	my $tamany_finestra = $_[2];		### LLargada del lloc de splicing
	my $posexo = $_[3];
	my $posreal = $_[4];
	my $tamany_fila = $_[5];
	my @puntuacio_nt;			### Vector on enregistrarem les puntuacions
	my $suma;				### Variable on enregistrem la puntuacio per cada nucleotid de que sigui un lloc de splicing
	my $x;					### Conte la posicio de la matriu de pesos
	my $n = 0;
	my $inici = $posreal - ($posexo-1);
	
	while ($n < ($tamany_fila - 1)) {	    ### Mentres podem fer finestres de llargada $tamany_finestra
		$suma = 0;
		$x = 0;    
		while ($x < $tamany_finestra) {			    ### Mentres no arribem al tamany de la finestra
			$suma = $suma + $matriu{$seq[$inici+$n+$x]}[$x];  ### Anem sumant les puntuacions 
			$x = $x + 1;
			}
		$puntuacio_nt[$n] = $suma;			    ### Enregistrem la puntuacio per aquell nucleotid
		$n = $n + 1;			
		}
	return(@puntuacio_nt);
}





#################### SOBRE_LLINDAR ################################ 
###Aquesta funcio ens selecciona els llocs de splicing acceptors i donadors que superen el llindar introduit
### Parametres d'entrada:
######################### Adreca del vector que conte a cada posicio la puntuacio d'una finestra d'una fila
######################### Valor del llindar a partir del qual considerem un lloc de splicing candidat com a lloc de splicing real
######################### Posicio del lloc de splicing en la que comenca o acaba a codificar per un exo
### Retorna: Matriu on tenim enregistrades les dades dels llocs de splicing que han superat el llindar
####################################################################


sub SOBRE_LLINDAR {
	my @v = @{$_[0]};  	### Vector que conte a cada posicio la puntuacio d'una finestra d'una fila
	my $llindar = $_[1];
	my $c = 0; 		### Comptador de les posicions del vector
	my @taula; 		### Matriu on enregistrarem les dades dels llocs de splicing que han superat el llindar
	my $f;
	### Numero a sumar al nucleotid on s'inicia (per acceptors) o acaba (per donadors) el lloc de splicing per tenir be l'inici de la regio codificant de l'exo
	
	while ($c < scalar(@v)) {				### Recorrem tot el vector
			if ($v[$c] >= $llindar) {		### Si aquell valor supera el llindar
				$taula [$f][0] = $c;		### En la primera columna enregistra el nucleotid que pot ser inici o fi d'un exo
				$taula [$f][1] = $v[$c];	### En la segona columna enregistra la puntuacio del lloc de splicing
				$f = $f +1;
				}
			$c = $c + 1;
			}
	return(@taula);
}


#################### MATRIU_ATG ################################### 
###Aquesta funcio ens busca els codons d'inici ATG, enregistra les posicions en una matriu. A la primera columna de la matriu tenim la posicio del primer nucleotid del codo ATG, i a la segona columna un zero. Aixo ho fem aixi per a que quan despres combinem les matrius, les tinguem totes en el mateix format. 
### Parametres d'entrada:
######################### Adreca del vector que conte una sequencia de nucleotids on a cada posicio tenim un nucleotid, i els ultims tres nucleotids formen un codo STOP (exceptuant la ultima fila)
######################### Pauta de lectura
######################### Fila de la matriu en la que ens trobem
### Retorna: Matriu on tenim enregistrades les dades dels codons ATG
####################################################################

sub MATRIU_ATG {

	my @fila = @{$_[0]};
	my $pauta = $_[1];
	my $num_fila = $_[2];		### Conte el numero de la fila a la que pertany la sequencia introduida
	my $inicial = 0; 		### Inici de la sequencia
	my $f = 0; 			### Comptador de les files de la matriu ATG
	my @matriu_ATG; 		### Matriu on posarem a la primera columna la posicio del codo ATG i a la segona un zero
	my $i = 0;			### Nucleotid pel que comencem a buscar
	
	
	if ($num_fila == 0) {		### Nomes en el cas de que ens entrin la primera fila, sumarem la pauta per ajustar la divisio dels codons
		$i = 0 + $pauta;
	} 							

### Busquem una sequencia de tres nucleotids formada per ATG en la sequencia donada

	while ($i<scalar(@fila)) {
		if ($fila[$i] =~ m/A/) { 
			if ($fila[$i+1] =~ m/T/){
				if ($fila[$i+2] =~ m/G/){
					$matriu_ATG[$f][0] = $i;
					$matriu_ATG[$f][1] = 0;
					$f = $f + 1;
					}
				}
			}
		$i = $i + 3;		
	}

	return(@matriu_ATG);
}



#################### MATRIU_STOP ################################### 
###Aquesta funcio ens mira si els tres ultims nucleotids de la sequencia introduida formen un codo STOP. Si es aixi, enregistra les posicions en una matriu. A la primera columna de la matriu tenim la posicio del tercer nucleotid del codo STOP, ja que aquest pot determinar la fi d'un exo. A la segona columna enregistrem un zero. Les dades les enregistrem en una matriu, encara que com a molt nomes pot trobar un codo stop per sequencia introduida. Aixo ho fem aixi per a que quan despres combinem les matrius, les tinguem totes en el mateix format. 
### Parametres d'entrada:
######################### Adreca del vector que conte una sequencia de nucleotids on a cada posicio tenim un nucleotid, i els ultims tres nucleotids formen un codo STOP (exceptuant la ultima fila)
### Retorna: Matriu on tenim enregistrades les dades dels codons STOP
####################################################################


sub MATRIU_STOP {
	my @s = @{$_[0]};		
	my $ultim_codo;			### En aquesta variable crearem l'ultim codo
	my @taula_STOP;			### Matriu on enregistrarem la taula amb les posicions dels codons STOP
	
	$ultim_codo = $s[(scalar(@s)-3)].$s[(scalar(@s)-2)].$s[(scalar(@s)-1)];		### Concatanem els tres darrers nucleotids per tal de crear un codo
	
	if ($ultim_codo =~ m/(TAA|TAG|TGA)/) {						### Si aquest ultim codo correspon a un codo STOP, enregistrem la posicio de l'ultim nucleotid del codo en la primera columna de la matriu, i en la segona columna hi enregistrem un zero
		$taula_STOP[0][0] = (scalar(@s)-1);
		$taula_STOP[0][1] = 0;
		}
	return(@taula_STOP);
}		



#################### COMBINACIONS ################################## 
###Aquesta funcio fa les combinacions compatibles donades dues matrius amb llocs d'inici (codo ATG o una sequencia donadora de lloc de splicing) i de fi d'exons (codo STOP o una sequencia acceptora de lloc de splicing). 
### Parametres d'entrada:
######################### Adreca de la matriu d'inici d'exo
######################### Adreca de la matriu de final d'exo 
######################### Paraula que defineix el tipus d'exo que creem a l'efectuar les combinacions (first, internal, terminal, single)
### Retorna: Matriu on tenim enregistrades les dades de les combinacions realitzades
####################################################################

sub COMBINACIONS {
	my @taula_inici = @{$_[0]};		### Matriu d'inici d'exo
	my @taula_fi = @{$_[1]};		### Matriu de final d'exo 
	my $tipus = $_[2];			### Paraula que defineix el tipus d'exo que creem a l'efectuar les combinacions
	my @taula_combinacions;			### Matriu on tenim enregistrades les dades de les combinacions realitzades
	my $fa = 0; 				### Comptador de les files de la taula d'acceptors
	my $fd; 				### Comptador de les files de la taula de donadors
	my $cf; 				### Comptador de les files de la taula de combinacions
	
	while ( $fa < scalar(@taula_inici)) {						### Recorrem tota la taula d'inici
		$fd = 0;	
		while($fd < scalar(@taula_fi)) {					### Recorrem tota la taula fi
			if ($taula_inici [$fa][0] < ($taula_fi [$fd][0]-9)) {		### En el cas de que el nucleotid d'inici sigui menor al de fi fer:
				$taula_combinacions[$cf][0] = $taula_inici[$fa][0];	### A la primera columna enregistrem el nucleotid d'inici de l'exo
				$taula_combinacions[$cf][1] = $taula_fi[$fd][0];	### A la segona columna enregistrem el nucleotid de fi de l'exo
                                $taula_combinacions[$cf][2] = $taula_inici[$fa][1];    	### A la tercera columna enregistrem la puntuacio del lloc acceptor de splicing (donat que el codo ATG no te cap puntuacio pero li assignem un valor de zero) 
				$taula_combinacions[$cf][3] = $taula_fi[$fd][1];	### A la quarta columna enregistrem la puntuacio del lloc donador de splicing (donat que el codo STOP no te cap puntuacio pero li assignem un valor de zero).
				$taula_combinacions[$cf][4] = $tipus;			### A la cinquena columna tenim el tipus d'exo
				$cf = $cf +1;
			}
			$fd = $fd +1;
		}	
		$fa = $fa +1;
		}
	return(@taula_combinacions);
}


#################### BIAIX_CODIFICANT ################################# 
###Aquesta funcio calcula el biaix codificant per tots els possibles exons que hem anotat a la matriu de combinacions
### Parametres d'entrada:
######################### Adreca de la matriu de combinacions
######################### Adreca del Hash en el que tenim la taula amb les frequencies d'us dels codons
######################### Adreca d'una sequencia de nucleotids
######################### Pauta de lectura. Aquesta no es la pauta de lectura de l'exo, sino que es la pauta amb la que hem generat la divisio d'exons al principi
### Retorna: Matriu amb les puntuacions i caracteristiques de cada exo
######################################################################

sub BIAIX_CODIFICANT {

my @combinacions = @{$_[0]};		### Matriu de combinacions		
my %hash = %{$_[1]};			### Hash en el que tenim la taula amb les frequencies d'us dels codons
my @s = @{$_[2]}; 			### Sequencia original d'una fila
my $pauta = $_[3];			### Pauta de lectura
my $num_fila = $_[4];

my $f = 0; 				### Comptador de les files de taula de combinacions
my $inici; 				### Nucleotid a partir del que comencem a calcular el biax codificant 
my $codo; 				### Enregistrarem cada cop 3 nucleotids concatenats
my $n; 					### Comptador del vector seq_interna
my @seq_interna; 			### Vector de la sequencia on posem a cada posicio 3 nucleotids
my $suma;	 			### Enregistrem la suma de les puntuacions segons el biaix
my @taula_puntuacio_files; 		### Enregistrem a taula la puntuacio de cada fila
my $resultat; 				### Aqui tindrem la suma de les puntuacions d'splicing i la del biaix codificant
my $pauta_real;

if ($num_fila != 0) {
	$pauta = 0;
}

while ($f < scalar(@combinacions)) {			### Per tots els intervals de la sequencia (exons)
	
	$inici = $combinacions[$f][0];		### Li assignarem el valor del lloc donador
	$suma = 0.0;
	
### La pauta real de l'exo sera:
	### 0: si la primera base de l'exo correspon a la primera base del codo al qual pertany (segons la divisio incial de la sequencia en codons)
	### 1: si la primera base de l'exo es la segona base del seu codo (segons la divisio incial de la sequencia en codons)
	### 2: si la priemra base de l'exo es la tercera base del seu codo (segons la divisio incial de la sequencia en codons)


	$pauta_real = 3; 				### Inicialment a la pauta de l'exo li assignem un valor de 3
	
	while(($inici - $pauta) % 3 != 0) {		### Si el nucleotid d'inici menys la pauta per la que vam fer la divisio de codons no es divisible per 3
		$inici = $inici + 1;			### Aleshores li sumem 1 al nucleotid d'inici
		$pauta_real = $pauta_real - 1;		### I a la pauta real li restem 1
	} 
	
	if ($pauta_real == 3) {				### Si la pauta real no ha estat modificada
		$pauta_real = 0;			### Vol dir que estem en la pauta 0
	}
	
	while ($inici <= (($combinacions[$f][1])-2)) {	### Mentres poguem fer codons comencant pel nucleotid d'inici 
		$codo = $s[$inici].$s[$inici+1].$s[$inici+2];
		$suma = $suma + log($hash{$codo}) - log(1/64); ### Com serien numeros molt petits, treballem amb el logaritme.
		$inici= $inici + 3;
		}

### Ara tenim el vector @seq_interna on a cada posicio te un codo. Ara buscarem cada codo com a etiqueta del hash de biaix codificant i anirem sumant la rao de versemblanca de la proporcio de cada codo

	$resultat = $combinacions[$f][2] + $combinacions[$f][3] + $suma;  ### Calculem la suma de les puntuacions d'splicing i la del biaix codificant
	$taula_puntuacio_files[$f][0] = $combinacions[$f][0];		  ### A la primera columna posem el nucleotid d'inici de l'exo
	$taula_puntuacio_files[$f][1] = $combinacions[$f][1];		  ### A la segona columna posem el nucleotid final de l'exo
	$taula_puntuacio_files[$f][2] = $resultat;			  ### A la tercera columna posem la puntuacio
	$taula_puntuacio_files[$f][3] = $pauta_real;			  ### A la quarta columna anotem la pauta de l'exo
	$taula_puntuacio_files[$f][4] = $combinacions[$f][4];		  ### A la cinquena columna anotem el tipus d'exo
	$taula_puntuacio_files[$f][5] = $suma;
	$f = $f + 1;
	}


	
return (@taula_puntuacio_files);
}




#################### POSICIONS_REALS ################################# 
###Aquesta funcio tradueix les posicions dels nucleotids de la matriu de resultats, en les posicions de la sequencia original. Aquesta traduccio cal fer-la, ja que les posicions dels nucleotids en la matriu les tenim respecte a la fila que hem creat, i per tant no equival a la posicio de la sequencia inicial.
### Parametres d'entrada:
######################### Adreca de la matriu de resultats
######################### Adreca de la matriu on a cada fila te una sequencia que acaba en codo STOP
######################### Numero de la fila que volem trobar la posicio real d'inici i fi
######################### Strand al que pertanyen aquests exons
### Retorna: Matriu amb les puntuacions i caracteristiques de cada exo, amb les posicions corregides
######################################################################


sub POSICIONS_REALS {

my @taula_puntuacio_files = @{$_[0]}; 		### Conte la matriu de resultats
my @matriu_seq = @{$_[1]}; 			### Matriu que a cada fila te una sequencia que acaba en codo STOP. Ho farem per cada pauta 
my $fila = $_[2]; 				### Numero de la fila que volem trobar la posicio real d'inici i fi
my $strand = $_[3];

my $c = 0; 	
my @taula_posicio_real; 			### Matriu de resultats corregida
my $p = 0; 					### Valor que hem de sumar per corregir les posicions			
my $j = 0;


	while ($c < $fila) {								### Mentres no arribem a la darrera fila
		$p = $p + scalar(@{$matriu_seq[$c]});					### Anirem sumant tots els nucleotids d'aquella fila
		$c = $c + 1;
		}
	
	while ($j < scalar(@taula_puntuacio_files)) {					### Per totes les files de la matriu de resultats fer:
		$taula_posicio_real[$j][0] = $taula_puntuacio_files[$j][0] + $p; 	### El nucleotid d'inici amb la posicio corregida
		$taula_posicio_real[$j][1] = $taula_puntuacio_files[$j][1] + $p; 	### El nucleotid final amb la posicio corregida
		$taula_posicio_real[$j][2] = $taula_puntuacio_files[$j][2]; 		### La tercera columna tindra el resultat de la puntuacio 
		$taula_posicio_real[$j][3] = $taula_puntuacio_files[$j][3]; 		### A la quarta columna anotem la pauta de l'exo
		$taula_posicio_real[$j][4] = $taula_puntuacio_files[$j][4]; 		### A la cinquena columna anotem el tipus d'exo
		$taula_posicio_real[$j][5] = $strand;				### Anotem l'strand
		$taula_posicio_real[$j][6] = $taula_puntuacio_files[$j][5];
		
		
		$j = $j + 1;
		}

return (@taula_posicio_real);
}



#################### RESULTATS_JUNTS ################################# 
###Aquesta funcio ajunta tres matrius
### Parametres d'entrada:
######################### Adreca de la primera matriu
######################### Adreca de la segona matriu
######################### Adreca de la tercera matriu
### Retorna: Matriu amb les tres matrius ajuntades
######################################################################

sub RESULTATS_JUNTS {

	my @R0 = @{$_[0]}; 	### Primera matriu
	my @R1 = @{$_[1]}; 	### Segona matriu
	my @R2 = @{$_[2]}; 	### Tercera matriu

	my @resultats_junts; 	### Matriu on les ajuntem

	my $c; 			### Comptador que compta les files de les matrius
	my $s = 0; 		### Comptador que compta les files de la matriu on ajuntem
	

### Hem hagut de fer els tres whiles's d'acontinuacio degut a que la longitud (en files) de cada matriu pot no ser la mateixa

	$c = 0;
	while ($c < scalar(@R0)) { 			### Per totes les files de la primera matriu 
				
		$resultats_junts[$s] = $R0[$c]; 	### Les copiem a una nova matriu
						
		$c = $c + 1;
		$s = $s + 1;
		}

	$c = 0;
	while ($c < scalar(@R1)) {			### Per totes les files de la segona matriu 
		$resultats_junts[$s] = $R1[$c];		### Les copiem a continuacio de la matriu que estem creant
			
		$c = $c + 1;
		$s = $s + 1;
		}

	$c = 0;
	while ($c < scalar(@R2)) {			### Per totes les files de la tercera matriu 
		$resultats_junts[$s] = $R2[$c];		### Les copiem a continuacio de la matriu que estem creant
		
		$c = $c + 1;
		$s = $s + 1;
		}


return(@resultats_junts);

}


#################### EXONS_BONS ###################################
### Aquesta funcio ens fa un cribatge de tots els exons predits, seleccionant nomes aquells que superin el llindar fixat per l'usuari
### Parametres d'entrada:
######################### Adreca de la matriu de resultats
#######a################## LLindar fixat per l'usuari
### Retorna: Matriu de resultats amb nomes els exons seleccionats
####################################################################

sub EXONS_BONS {

my @matriu = @{$_[0]};			### Matriu de resultats
my $llindar = $_[1];			### LLindar fixat per l'usuari
my $c = 0; 				### Comptador de files de la matriu de resultats
my $fila = 0; 				### Comptador de files de la matriu amb els resultats seleccionats
my @matriu_seleccionada;		### Matriu amb els resultats seleccionats

 
	while ($c<scalar(@matriu)) {				     ### Per totes les files de la matriu
		if  ($matriu[$c][2] >= $llindar) {		     ### Si la puntuacio de l'exo supera el llindar fer:
			$matriu_seleccionada[$fila] = $matriu[$c];   ### Copiem la fila en la nova matriu
			$fila = $fila + 1; 
			}
		$c = $c + 1;
		}
return (@matriu_seleccionada);
}


#################### POSICIONS_REALS_STRANDNEG ################################# 
###Aquesta funcio calcula les posicions reals dels nucleotids que formen l'exo en l'strand negatiu. Entenem com a posicions reals, les posicions en la sequencia introduida per l'usuari.
### Parametres d'entrada:
######################### Adreca de la matriu de resultats
######################### Longitud de la sequencia inicial
### Retorna: Matriu on tenim enregistrades les dades dels exons amb les posicions reals dels nucleotids de l'strand negatiu
################################################################################

sub POSICIONS_REALS_STRANDNEG {

my @resultat_negatiu =  @{$_[0]}; 		### Conte la matriu de resultats de l'strand -
my $scalar = $_[1]; 				### Scalar de la sequencia original

my @resultat_neg_real;				### Matriu on enregistrarem els resultats amb les posicions correctes
my $c = 0; 					### Comptador de les files de la nova matriu


	while ($c < scalar(@resultat_negatiu)) {					### Mentres hi hagi files a la matriu, fer:

		$resultat_neg_real[$c][1] = ($scalar - 1) - $resultat_negatiu[$c][0];   ### Inici tindra el valor a strand + i sera la posicio 1 ja que es el mes gran
		$resultat_neg_real[$c][0] = ($scalar - 1) - $resultat_negatiu[$c][1];   ### Fi tindra el valor a strand + i sera la posicio 0 ja que es el mes petit
		$resultat_neg_real[$c][2] = $resultat_negatiu[$c][2];		        ### La tercera columna tindra el resultat de la puntuacio que ja tenia l'strand
		$resultat_neg_real[$c][3] = $resultat_negatiu[$c][3]; 			### A la quarta columna anotem la pauta de l'exo
		$resultat_neg_real[$c][4] = $resultat_negatiu[$c][4]; 			### A la cinquena columna anotem el tipus d'exo
		$resultat_neg_real[$c][5] = $resultat_negatiu[$c][5];			### Strand negatiu
	

	$c= $c + 1;
		}

return (@resultat_neg_real);
}


#################### UNIR_RESULTATS ################################# 
####Aquesta funcio ajunta dues matrius
### Parametres d'entrada:
######################### Adreca de la primera matriu
######################### Adreca de la segona matriu
### Retorna: Matriu amb les dues matrius ajuntades
#####################################################################

sub UNIR_RESULTATS {

my @pos = @{$_[0]};  		### Primera matriu
my @neg = @{$_[1]};		### Segona matriu
my @unit; 			### Matriu amb els resultats units
my $fila = 0; 			### Comptador de les files de cada matriu 
my $c = 0; 			### Comptador de les files de la nova matriu @unit

while ($fila < scalar(@pos)) {	### Mentres hi hagi files en la matriu, fer:
	$unit[$c] = $pos[$fila];### Anem copiant les files en la nova matriu
	$fila = $fila + 1;
	$c = $c + 1;
	}

$fila = 0; 

while ($fila < scalar(@neg)) {  ### Mentres hi hagi files a la matriu, fer:
	$unit[$c] = $neg[$fila];### Anem copiant les files a continuacio de la matriu que estem creant
	$fila = $fila + 1;
	$c = $c + 1;
	}

return(@unit);
}


#################### FORMAT_GFF ################################### 
### Aquesta funcio crea un fitxer de sortida en format GFF a partir de les dades dels exons predits, enmagatzemades en una matriu
### Parametres d'entrada:
######################### Identificador de la sequencia
######################### Adreca de la matriu amb les dades dels exons predits
######################### Nom del fitxer de sortida que crearem
### Retorna: No retorna res. Aquesta funcio simplement l'hem feta per agrupar un conjunt d'accions.
####################################################################

sub FORMAT_GFF {

my $i = $_[0]; 			### Identificador
my @matriu = @{$_[1]}; 		### Matriu de resultats
   
my $c = 0; 			### Comptador de les files de la matriu


### Escriurem en el fitxer de sortida la seguent informacio:
print "#Exons predits amb el CERCAEXONS per la sequencia $i\n\n";

while ($c < scalar(@matriu)) {			### Mentres hi hagi files a la matriu, fer:
  	print "$i\t"; 					### La primera columna tindra l'identitat de la sequencia 
	print "cercaexons\t"; 				### La segona columna tindra el nom del nostre programa
	print "$matriu[$c][4]\t"; 			### La tercera columna especificara el tipus d'exo : first, internal, terminal o single
	print $matriu[$c][0]+1,"\t";			### La quarta columna tindra la posicio de l'inici de l'exo en la sequencia
	print $matriu[$c][1]+1, "\t"; 			### La cinquena columna tindra la posicio de la fi de l'exo en la sequencia
	printf("%.2f \t", $matriu[$c][2]);		### La sisena columna tindra la puntuacio de l'exo
	print "$matriu[$c][5]\t"; 			### La setena columna indicara l'estrand en que es llegeix l'exo 
	print "$matriu[$c][3]\n"; 			### La vuitena columna tindra indicara la pauta de lectura

	$c = $c + 1;
	}

close(GFF);
}

################### FI DEL PROGRAMA ##############################################################

