Materials i mètodes



AUTOMATITZACIÓ

Per tal d'automatitzar els programes que havíem de fer servir, hem creat uns sèrie de programes en llenguatge bash i en llenguatge perl per tal de facilitar el processament d'un gran volum de dades. Concretament hem creat 5 programes bash i 2 programes perl:

⁃ tblastnhs.bash
⁃ classificador.pl
⁃ exoneratehs.bash
⁃ exchange.pl
⁃ genewisehs.bash
⁃ tcoffee.bash
⁃ SECIS.bashb

Per tal que algun usuari desitgi fer-los servir hem escrit un tutorial bàsic:

En primer lloc, hem de crear una carpeta al directori dels programes anomenada seqaahs, a on s'han de col·locar totes les seqüències que s'utilitzaran com a query en un arxiu diferent en format fasta (ex: DI2.fa). Per altra banda, hem de crear una variable que contingui el path fins el nostre genoma ($panda = /cursos/BI/genomes/project_2013/Ailuropoda_melanoleuca/genome.fa).Per últim, cal tenir instal·lats els programes blast, exonerate, genewise i t_coffee.

1. tblastnhs.bash: és un programa que produeix alineaments locals de seqüències tipus tBLASTn de forma automàtica. Permet comparar un conjunt de seqüències proteiques dins de la carpeta seqaahs contra la base de dades de nucelòtids previament definida (el genoma d'Ailuropoda melanoleuca) al directori cursos/BI/genomes/project_2013/Ailuropoda_melanoleuca/genome.fa.

Requereix:
· Una carpeta al directori anomenada ouput.blast, on s'enviaran els fitxers de sortida dels blasts amb el nom “query.hs.txt”.
· Una segona carpeta anomenada blast.selec: és el directori on es generarà un arxiu per a cada hit trobat amb una e-value significativa (e-value < 0,001) amb el nom "query_hit.txt".
· El programa perl “classificador.pl”.

Per visualitzar i valorar el resultat obtingut es crearà un arxiu anomenat "resultats_blast.txt" on es veuran només els hits que siguin significatius.

Codi BLAST: 

#!/bin/bash

#$ -o tblastn.stdout
#$ -e tblastn.stderr
#$ -q llicen.q
#$ -N blastJob
#$ -cwd

source ~/.bashrc

### Creem una variable que contingui el path fins el genoma del panda

panda="/cursos/BI/genomes/project_2013/Ailuropoda_melanoleuca"

### Això és per a detectar quines proteines hi ha a la carpeta seqaahs

for proteines_fa in ./seqaahs/*; do {
basename_fa=`basename $proteines_fa`
proteina=${basename_fa%.fa}
proteines=$proteines" $proteina"
} done

### Programa: fem dos bucles for per tal de fer les combinaciones entre proteïnes (query) i el genoma del panda


	for protein in $proteines ; do {
		
			patroalfa=`cat ./seqaahs/$protein.fa | egrep ^">"`
			patrobeta=${patroalfa#*>}
			echo $patrobeta > temppatro
			patro=`cut -f 1 -d " " temppatro`
			rm temppatro
			
			echo ">${protein}:" >> resultats_blastsh.txt
			
			blastall -p tblastn -i ./seqaahs/$protein.fa -d /cursos/BI/genomes/project_2013/Ailuropoda_melanoleuca/genome.fa -m9 -o output.blast/$protein.txt 
			
			x=`egrep $patro ./output.blast/$protein.txt | egrep e-`
			if [ "$x" != "" ]; then
				egrep $patro output.blast/$protein.txt | egrep e- > arxiu.temp
				mkdir ./blast.selec/$protein 2> error.temp
				./classifica.pl $protein < arxiu.temp

			fi
		
		
			echo ${protein} fet
		} done
	} done

rm ./*.temp

Codi classificador:

#!/usr/bin/perl

use strict;
my $contador=0;

# Amb aquest programa classifica els resultats del tblastn en carpetes amb el nom de cada proteïna analitzada, col·locant dintre els hits obtinguts

while ()
{
	open (fitxersortida, ">./blast.selec/$ARGV[0]/hit$contador.txt");
	print fitxersortida $_;
	close fitxersortida;
	$contador++;
}

2. exoneratehs.bash: Aquest programa executarà l'exonerate per als hits obtinguts amb el programa "tblastnhs.bash".

Requereix la creació de les següents carpetes en el mateix directori que el programa:
· “subseqshs” on es guardaran totes les regions extretes mitjançant fastasubseq.
· “exonerate_ouput” és el directori on s'enviaran els fitxers de sortida de l'exonerate en format gff.
· "cDNA". En aquest directori es guardaran els cDNAs de les proteïnes predites per l'exonerate (en cas que hagi estat capaç de predir-ne).
· "translations". es guardarà un arxiu amb les 6 possibles traduccions per a cada cDNA obtingut.
· El programa perl “exchange.pl”.

Per saber quins exonerates han estat capaços de predir gens es crea un arxiu anomenat "resultats_exonerate.txt" on només es mostren els exonerates que han predit gens.

Codi exoneratesh.bash:

!/bin/bash

rm resultats_exonerate.txt 2> error.temp

#$ -o tblastn.stdout
#$ -e tblastn.stderr
#$ -q llicen.q
#$ -N blastJob
#$ -cwd

source ~/.bashrc

### Això és per poder triar si volem exonerate exhaustiu o no

if [ "$1" = "-exhaustive" ]; then
exhaustive="--exhaustive yes"
fi

## Per cada hit de cada proteïna dintre de la carpeta blast.selec:

for proteins in blast.selec/* ; do {
	for hit in $proteins/*; do {

		proteins=`basename $proteins`
        	hit=`basename $hit`
		hit=${hit%.txt}

		echo ${proteins}_${hit}:

		### Això és per a detectar l'identificador que apareix al fitxer fasta a la primera línea

		patro=`cut -f 2 blast.selec/$proteins/$hit.txt` 

		### Per emagatzemar la regió on apareix el hit

		fastafetch $panda panda.index $patro > $patro.fa

		### Aquí definim l'inici i el final del hit en la regió

		inici=`cat ./blast.selec/$proteins/$hit.txt | cut -f 9`
		final=`cat ./blast.selec/$proteins/$hit.txt | cut -f 10`

		if [ "$inici" -gt "$final" ]; then
			inici="$final"
		fi
		
		start=$(( $inici - 10000 ))
		length="30000"

		### Amb al comanda fastasubseq extraiem una regió del genoma del panda		

		fastasubseq $patro.fa $start $length > ./subseqhs/${proteins}.${hit} 2> error_msg.temp
x
		### Programa exonerate per a la predicció de gens
		
		echo "	Iniciant exonerate:"

		./exchange.pl < ./seqaahs/$proteins.fa > ./seqaahsx/$proteins.fa

		exonerate -m p2g --showtargetgff -q ./seqaahsx/$proteins.fa -t ./subseqhs/${proteins}.${hit} $exhaustive  > ./exonerate_output.hs/${proteins}_${hit}.exonerate 2> error_de_substitucio.temp
				
		gene=`grep gene ./exonerate_output.hs/${proteins}_${hit}.exonerate`
		if [ "$gene" != "" ]; then
			echo ${proteins}_${hit} : >> exonerate_output.hs/resultats_exonerate.txt
		fi

		echo "	Exonerate acabat"	

		### Per extreure el cDNA i traduir-lo a proteïnes fem servir el programa fastaseqfromGFF.pl i la comanda fastatranslate				
		
		echo "	Preparant i executant per a t_coffee:"
				
		exon=`grep -w exon ./exonerate_output.hs/${proteins}_${hit}.exonerate`
		grep -w exon ./exonerate_output.hs/${proteins}_${hit}.exonerate > exonerate_output.hs/${proteins}.${hit}.gff
		if [ "$exon" != "" ]; then
			fastaseqfromGFF.pl ./subseqhs/${proteins}.${hit} exonerate_output.hs/${proteins}.${hit}.gff > ./cDNA.hs/${proteins}_${hit}.fa
			fastatranslate ./cDNA.hs/cDNA_${proteins}_${hit}.fa > ./translations.hs/aa_${proteins}_${hit}.fa
		fi
	
		echo ${proteins}_${hit} fet
		
		} done
	} done

rm *.temp

Codi exchange.pl:

#!/usr/bin/perl
use strict;

my @v;
my $string;
my $counter = 0;

# Amb aquest programa canviem les U de les proteïnes.fa per X per tal de porder realitzar el prorgama exoenrate

while ()
{
	if ($counter == 0)
	{
		print $_;
	}
	if ($counter > 0)
	{

		@v=split(//,$_);
	
		while (<@v>)
		{
			if ($_ eq "U")
			{		
				$string.= "X";
			}		
			else
			{
				$string.=$_;
			}
		}
	}
	$counter++;
	$string.="\n";
}

print $string;

3. genewisehs.bash: Aquest programa executarà el genewise per als hits obtinguts amb el programa "tblastn.bash".

Requereix:
· Una carpeta al mateix directori anomenada “genewise_ouput” on s'enviïn els arxius de sortida del genewise.

Codi genewisehs.bash:

#!/bin/bash

# Programa genewise per a la predicció de gens

for proteins in blast.selec/* ; do {
        for hit in $proteins/*; do {

		proteins=`basename $proteins`
                hit=`basename $hit`
                hit=${hit%.txt}

		genewise -pep -pretty -cdna -gff -both ./seqaahsx/$proteins.fa ./subseqhs/${proteins}.${hit} > ./genewise_output.hs/${proteins}_${hit}.txt		
s		
		echo ${proteins}_${hit} fet

		} done
	} done

4. tcoffee.bash: Aquest programa executarà un alineament amb T-Coffee amb els resultats de l'exonerate i el genewise amb la seva respectiva query.

Requereix la creació de les següents carpetes en el mateix directori que el programa:
· “lectures” on es trobin els marcs de lectura 1 de les proteïnes predites.
· “tcoffee.aln” i “tcoffee.html” on s'enviaran els fitxers de sortida .aln i .html respectivament.

1. Per saber quins tcoffee s'han dut a terme es crea un arxiu anomenat “tcoffee.resultats” on nomes és mostren els noms de les querys amb els que s'ha fet el tcoffee.

Codi tcoffee.bash:

!/bin/bash

# Programa t_coffee per l'alineament de proteïnes

for file in lectures/* ; do {
	name=`basename $file`
	name=${name%.fa}
	protein=${name%.hit*}

echo $name

echo $protein

	t_coffee seqaamm/$protein.mm.fa $file > ./tcoffee.resultats/tcoffee.$name.txt 2> error.temp
	mv $protein.mm.aln ./tcoffee.aln/tcoffee.$name.aln
	mv $protein.mm.html ./tcoffee.html/tcoffee.$name.html

echo $file done

} done

5. SECIS.sh: Aquest programa cercarà les possibles estructures SECIS que es troben en els subseqs de les diferents proteïnes.

Requereix:
· Creació d'una carpeta en el mateix directori anomenada “secis_ouput” on s'enviaran els fitxers de sortida del programa.

Codi SECIS.sh:

#!/bin/bash

#Programa per automatizar la busca d'elemets SECIS amb el programa SECIsearch;

echo "A la recerca d'elements SECIS";

export PATH=/cursos/BI/bin:$PATH;

for proteins in subseq/* ; do {
	for hit in $proteins/*; do {
		proteins=`basename $proteins`
                hit=`basename $hit`
		hit=${hit%.}
			                        
	SECISearch.pl ./subseq/${proteins}/${hit} > ./secis_output/secis_${hit}.txt

echo ${hit}.secis fet

	} done 
} done
Up