Construction hash intron

**shadow19c** · 10/03/2011, 17h12

Bonjour,
voila toujours des problemes avec le perl,
mon souci c'est que je dois construire un hash de taille d'intron a partir de fichier cds:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
 
# Colonne 1: scaffold
#         2: gene id
#         3: start
#         4: dir
#         5: origine gene model
#         6: proteinId
#         7: ex-nom
#         8: coordonnees genomiques
#         9: structure intron/exon
#        10: segments du CDS couvert par EST
C169-scaffold_1	C169v2-00001	3777	-	ORIGINAL JGI	55057	Genemark1.1_g	3777..3857,4046..4192,4443..4561,4940..5234,5406..5540,5734..5847,6009..6098,6421..6492
C169-scaffold_1	C169v2-00002	12682	-	ORIGINAL JGI	55058	Genemark1.2_g	12682..12691,13195..13445,13694..13711
C169-scaffold_1	C169v2-00003	18095	+	ORIGINAL JGI	31905	fgenesh1_kg.1_#_1_#_4092_1_CBOZ_CBPA	18095..18097,18280..18410,18690..18972
C169-scaffold_1	C169v2-00004	20452	+	ORIGINAL JGI	6968	gw1.1.615.1	20452..20496,20636..20726,20881..21046,21194..21382,21567..21735

Je dois recuperer les coordonnees et construire lhash pour chaque id (CV-****),
pour cela si j'ai une coordonne que voici: a..b,c..d,e..f , les "," correspondent aux introns pour cela donc je dois prendre le couple(b,c) et (d,e) et faire la soustraction du start et end de lintron et calculer par la suite la taille de l'intron.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
 
#!/usr/local/bin/perl
#use strict;
#use warnings;
 
open (DESC1, "coco_cds");
while(<DESC2>){
	next if /^\#/;
	chomp;
	my @exon=(split/,/,$coord);
	(my @intron)=($coord=~s/(/d+,/d+)/g);
	if ($strand eq '-'){
		reverse(@exon,@intron);
	}
	for(0<$i<întron){
		$exon=$exon[i];
		$intron=$intron[i];
		($s,$e)=(split/,/,$intron);
		$taille=abs($e-$s+1);
		push @cds,($exon);
		$aa=int(@cds/3);
		hash{$id}{$aa}=#taille je ne sais pas comment lui dire
	}
}

Voici mon code du debut mais j'arrive pas a regler le probleme du regex, mais aussi a lui dire que hash{$id}{$aa}=taille de l'intron et pour ca je dois faire pour chaque couple d'intron donc je devrais avoir pour cette exemple deux tailles dintron du couple (b,c) et (d,e).

Est ce que vous pouvez m'aider pour mon code .

Merci

**Jasmine80** · 11/03/2011, 09h08

bonjour,

je n'ai pas compris la fin, donc il faudra modifier mon code :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
 
#!/usr/local/bin/perl
 
use strict;
use warnings;
 
 
 
open (DESC1, "coco_cds");
 
while (my $line = <DESC1>){
 
	next if ($line =~ m/^\#/);
	chomp $line;
 
	my ( $id, $strand, $coord) = $line =~ m/([\w-]+)\t[\w-]+\t\d+\t([-+])\tORIGINAL JGI\t\d+\t[\#\w\.]+\t([,\d\.]+)/; 
	my @exon = split (/,/, $coord);
 
	my @intron = ($coord =~ m/(\d+,\d+)/g);
 
	if ($strand eq '-'){
		@exon = reverse(@exon);
		@intron = reverse(@intron);
	}
 
	my $cds;
	my %hash; 
 
# 3777..3857,4046..4192,4443..4561,4940..5234,5406..5540,5734..5847,6009..6098,6421..6492
#          324
 
	for my $i (0..$#intron){
 
		my ($s,$e)=(split/,/,$intron[$i]);
		my $taille=abs($e-$s+1);
 
		my ($s2,$e2)=(split/\.\./,$exon[$i]);
		my $taille2 =abs($e2-$s2+1);
 
		$cds += $taille2;
		my $aa=int($cds/3);
 
		$hash{$id}{$aa}=$taille;
 
	}
 
}

**shadow19c** · 11/03/2011, 17h01

Bonjour,
voila le but de mon long programme cest de chercher la taille des introns orthologues:
Pour cel a je pars de deux fichiers Cds Aster et CV(comme d'hab lol), ensuite je construis un hash, apres cela je lis mon fichier de reciprocal best blast hits qui va me permettre davoir le couple Cv-***, Aster_****.
Ensuite je dois construire un pipeline qui va me recuperer pour chaque CV et aster de couple me creer un fichier fasta donc je vais faire un alignement blast de deux sequences.
A partir de ce fichier d'aligenement je vais parcourir les deux sequences jusqu ce que j'obtiens le bon intron calcule lors du hash.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
 
#!/usr/local/bin/perl
 
use strict;
use warnings;
 
 
 
open (DESC1, "coco_cds");
 
while (my $line = <DESC1>){
 
	next if ($line =~ m/^\#/);
	chomp $line;
 
	my ( $id, $strand, $coord) = $line =~ m/([\w-]+)\t[\w-]+\t\d+\t([-+])\tORIGINAL JGI\t\d+\t[\#\w\.]+\t([,\d\.]+)/; 
	my @exon = split (/,/, $coord);
 
	my @intron = ($coord =~ m/(\d+,\d+)/g);
 
	if ($strand eq '-'){
		@exon = reverse(@exon);
		@intron = reverse(@intron);
	}
 
	my $cds;
	my %hash; 
 
# 3777..3857,4046..4192,4443..4561,4940..5234,5406..5540,5734..5847,6009..6098,6421..6492
#          324
 
	for my $i (0..$#intron){
 
		my ($s,$e)=(split/,/,$intron[$i]);
		my $taille=abs($e-$s+1);
 
		my ($s2,$e2)=(split/\.\./,$exon[$i]);
		my $taille2 =abs($e2-$s2+1);
 
		$cds += $taille2;
		my $aa=int($cds/3);
		 $hash{$id}{$aa}=$taille."\n";
 
 
	}
 
}
 
open (DESC2, "aster_cds");
 
while (my $line = <DESC2>){
 
	next if ($line =~ m/^\#/);
	chomp $line;
 
	my ( $id, $strand, $coord) = $line =~ m/([\w-]+)\t[\w-]+\t\d+\t([-+])\tORIGINAL JGI\t\d+\t[\#\w\.]+\t([,\d\.]+)/; 
	my @exon = split (/,/, $coord);
 
	my @intron = ($coord =~ m/(\d+,\d+)/g);
 
	if ($strand eq '-'){
		@exon = reverse(@exon);
		@intron = reverse(@intron);
	}
 
	my $cds;
	my %hash; 
 
# 3777..3857,4046..4192,4443..4561,4940..5234,5406..5540,5734..5847,6009..6098,6421..6492
#          324
 
	for my $i (0..$#intron){
 
		my ($s,$e)=(split/,/,$intron[$i]);
		my $taille=abs($e-$s+1);
 
		my ($s2,$e2)=(split/\.\./,$exon[$i]);
		my $taille2 =abs($e2-$s2+1);
 
		$cds += $taille2;
		my $aa=int($cds/3);
 
			$hash{$id}{$aa}=$taille."\n";
 
 
	}
 
}
 
open (DESC3, "R_B_B_H_A_C.txt");
while (my $ligne1 = <DESC3>){
 
	my ($id1,$id2)=(split/\t/,$ligne1);
 
	system("formatdb -o -pT -i DB.pep");
	system("fastacmd -d DB.pep -s $id1>i");
	system("fastacmd -d DB.pep -s $id2>j");
	system("bl2seq -i i -j j -p blastp -o out.bl -e 1.e-5");
 
}

Voci mon code du debut de pipeline, mon soucis cest que j'arrive pas lui dire de faire le blast2seq que pour le premier couple ensuite je vais faire la recherche de parsing puis des que cest termine de recommence.

Merci

**Jasmine80** · 14/03/2011, 07h47

pourquoi avoir ajouter un \n?

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

$hash{$id}{$aa}=$taille."\n";

pour ton blast tu peux utiliser le module Bio::Tools::Run::StandAloneBlast

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
	my @params = ('program' => 'blastn');
	my $factory = Bio::Tools::Run::StandAloneBlast->new(@params);
	$factory->outfile($out_file);
	$factory->executable("bl2seq", "C:/BLAST/bin/bl2seq.exe");
	my $input = Bio::Seq->new(-id=>$id1 ,-seq=> $seq1);
	my $input2 = Bio::Seq->new(-id=>$id2 ,-seq=> $seq2);
	my $bl2seq_report = $factory->bl2seq($input, $input2);

**shadow19c** · 14/03/2011, 10h39

Bonjour,
le "\n" cetait juste pour moi quand je voulais verfier, mais je l'ai enleve, par contre il y a une erreur a cette ligne:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
my ( $id, $strand, $coord) = $line =~ m/([\w-]+)\t[\w-]+\t\d+\t([-+])\tORIGINAL JGI\t\d+\t[\#\w\.]+\t([,\d\.]+)/;

Car le id doit correspondre a Aster-** ou CV-** et non a scalfold

J'ai corrige en faisant

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
my ( $id, $strand, $coord) = $line =~ m/[\w-]+\t([\w-]+)\t\d+\t([-+])\tORIGINAL JGI\t\d+\t[\#\w\.]+\t([,\d\.]+)/;

Merci, en tout cas , par contre je ne peux utiliser Bioseq car ce n'est pas installe sur les machnes et je nai pas les droits donc je lance des commandes systemes , mais ca marche comme ca j'apprend perl en meme temps.

**Jasmine80** · 14/03/2011, 10h55

oui, au lieu de cette regexp, tu peux utiliser split sur \t comme l'avait fait remarquer Norore précédemment, ça a l'avantage d'être simple et général

Construction hash intron

Bioinformatique Perl

Vue hybride

Discussions similaires

Partager

Partager