view catchsequence/catchsequence.pl @ 0:c30eb2050ef5 draft default tip

Uploaded
author mgarnier
date Thu, 26 Aug 2021 13:41:33 +0000
parents
children
line wrap: on
line source

#!/usr/bin/perl

use strict;
use warnings;

#INPUTS_
# my $Result_RES = $ARGV[0];
my $sequences = $ARGV[0];

#OUTPUT_
#my $output = $ARGV[1];

#my @list_seq = split(/,/,$sequences);
my @list_seq = @ARGV;

##########################################################################################

#open (OUT, ">$output");
print "SEQUENCE\tRESISTANCE GENES\tPLASMIDS\tVIRULENCE GENES\tMLST NUMBER\tGENE(ALLELE)\n";

foreach my $sequence (@list_seq) { 
	my $Result_RES = `abricate --db resfinder $sequence > $sequence.RES.txt`; #appel système de la commande abricate avec la BDD ResFinder
	my $Result_PLA = `abricate --db plasmidfinder $sequence > $sequence.PLA.txt`; #appel système de la commande abricate avec la BDD PlasmidFinder
	my $Result_VIR = `abricate --db vfdb $sequence > $sequence.VIR.txt`;
	my $Result_MLST = `mlst $sequence > $sequence.MLST.txt`;

	open (RES, "$sequence.RES.txt");
	print "$sequence\t";

	while (<RES>) {

    	chomp();
    		if ($_ !~ m/^#/) {
    			my @infos = split(/\t/,$_);
    			my $geneRes = $infos[5]; # resistance gene name (ancienne valeur $infos[4])
    			my $identity = $infos[10]; # identity % (ancienne valeur $infos[9])
  
    			if ($identity > 90.000) {
        			print "$geneRes;";
    			}
    		}

	}

	close (RES);
	print "\t";


	open (PLA, "$sequence.PLA.txt") or die "could not open $!";

	while (<PLA>) {
    	chomp();
    		if ($_ !~ m/^#/) {
        		my @infos = split(/\t/,$_);
        		my $plasmid = $infos[5]; # plasmid name
        		my $identity = $infos[10]; # identity %
 
        		if ($identity > 90.000) {
                		print"$plasmid;";
        		}
    		}

	}
	close (PLA);
	print "\t";

	open (VIR, "$sequence.VIR.txt") or die "could not open $!";

	while (<VIR>) {
    	chomp();
    	
		if ($_ !~ m/^#/) {
        		my @infos = split(/\t/,$_);
        		my $geneVir = $infos[5]; # virulence gene name
        		my $identity = $infos[10]; # identity %
 
        		if ($identity > 80.000) {
                		print "$geneVir;";
        		}
    		}

	}
	close (VIR);
	print "\t";


	open (MLST, "$sequence.MLST.txt") or die "could not open $!";

        while (<MLST>) {
        chomp();

        	my @infos = split(/\t/,$_);
                my $numMLST = $infos[2];
            	print "$numMLST\t";
   		
		for (my $i=3; $i <= $#infos; $i++){
			print "$infos[$i];";
		}

        }
        close (MLST);
        print "\n";
}

#close (OUT);

unlink glob ('*.VIR.txt');
unlink glob ('*.PLA.txt');
unlink glob ('*.RES.txt');
unlink glob ('*.MLST.txt');