Discussion :

[shadow19c]

Bonjour, voila j'ai un petit projet d'aligenement de sequences en gros et comme je suis nouveau en python, et j'ai pas appris le perl je voulais savoir si il y en a prmis vous qui font du python et qui peuvent m'aider le long de mon periple snif snif


merci

[wiztricks]

Salut,
Le forum est généralement hanté par de bonnes volontés.
Mais, il faudrait être un peu plus précis sur les questions.
- W

[shadow19c]

voila je dois faire de la statistique de intron je m'explique:
J'ai un fcihier venant de blast

Code :

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

1
2
3
4
5
6
7
8
 
C169v2-04488	Aster-06787	187	262	561	636       
C169v2-04487	Aster-06725	92	235	216	357
C169v2-04480	Aster-06724	371	431	4	64
C169v2-04481	Aster-06792	2	447	41	488	
C169v2-03069	Aster-04749	21	283	2	265
C169v2-03068	Aster-01309	8	222	41	248
C169v2-03061	Aster-06908	100	396	399	708

je dois pour chaque couple CV-**** et Aster-****
aller transformer les coordonnees genomiques en proteiques:
1er fichier des CV-***

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
 
# 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
C169-scaffold_1	C169v2-00005	21893	-	ORIGINAL JGI	34662	fgenesh1_pm.1_#_3	21893..22051,22201..22293,22590..22747,22864..23102,23192..23364
C169-scaffold_8	C169v2-04488	528365	+	ORIGINAL JGI	63423	Genemark1.4206_g	528365..528435,528774..528938,529294..529389,529779..530001,530348..530487,530870..530990
C169-scaffold_8	C169v2-04487	523752	+	ORIGINAL JGI	42018	fgenesh1_pg.8_#_66	523752..523949,524383..524517,524819..524948,525371..525456,525763..525830,526605..526740,527270..527305
C169-scaffold_5	C169v2-03069	1412897	+	ORIGINAL JGI	41001	fgenesh1_pg.5_#_204	1412897..1412909,1413103..1413155,1413409..1413531,1413634..1413748,1414228..1414427,1414721..1414891,1415181..1415381,1415548..1415679

2 eme fichier ASTER-***

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
 
# 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
scaffold_00001	Aster-00001	17594	-	ORIGINAL JGI	14863	gw1.00001.44.1	17594..17852,17997..18023,18357..18572,18689..18759,18842..18967
scaffold_00001	Aster-00002	19928	-	ORIGINAL JGI	21312	e_gw1.00001.25.1	19928..20186,20453..20616,20744..20792,20913..21081,21507..21714
scaffold_00001	Aster-00003	31765	+	ORIGINAL JGI	14113	gw1.00001.37.1	31765..31855,32261..32514,32748..32831,32917..32998,33154..33362,33646..33753
scaffold_00001	Aster-00004	44240	-	ORIGINAL JGI	36055	fgenesh1_pg.00001_#_3	44240..44511,44793..44943,45050..45133,45256..45318,45686..45721,45782..45904,46025..46342,46459..46642,46822..46871
scaffold_00001	Aster-00005	48441	+	ORIGINAL JGI	21325	e_gw1.00001.70.1	48441..48583,48791..48920
scaffold_00101	Aster-04749	530046	+	ORIGINAL JGI	29654	e_gw1.00101.62.1	530046..530292,530743..530942,531299..531488,531763..532046
scaffold_00145	Aster-06787	986197	-	ORIGINAL JGI	41618	fgenesh1_pg.00145_#_90	986197..986425,986528..986722,986990..987130,988711..988816,990091..990181,990354..990470,990778..990805,991167..991259,991536..991663,991869..991991,992198..992278,992360..992594,993369..993481,993632..993868

Je dois par exemple pour le premier couple C169v2-04488 Aster-06787 aller verfier dans leur fcihier respectifs si ils existent ensuite verfier le brin si cest + cest bon si cest - renverser les coordonnes car lecture a lenvers.

Je mexplique mieux:

prenons un exemple C169v2-03069 Aster-04749

on a des match pour C169v2-03069 dans le premier fichier et Aster-04749 dans le deuxieme.

Code :

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

1
2
3
4
5
6
7
8
 
with open('coco_cds.txt') as fichier:
	for line in fichier:
		ch= line.split('\t')[1:2]
		if ch in lsite:
			brin=line.split('\t')[3:4]
			if '+'in brin:	
				print "ok"

meme code pour lautre.

Je dois par contre transformer les coordonnes genomiques en coordonnes proteiques :

Code :

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

1412897..1412909,1413103..1413155,1413409..1413531,1413634..1413748,1414228..1414427,1414721..1414891,1415181..1415381,1415548..1415679

cest a dire 1412897 ceci est le premier donc je le met a 1 ensuite je fais (1412909-1412897)/3 donc jai 4 donc le premier exon est de 1-4 ensuite il ya une virgule donc cest lintron donc je peux dire que jai un intron en 4 ainsi de suite

pour tout les coordonnes.

A la fin je dois faire un print pour dire C169v2-03069 3 introns Aster-04749 4 introns
on ma dit de mettre pour chaque couple dans un dictionnaire.

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
 
#! /usr/bin/python
 
 
#coding=utf-8
coco=[]
aster=[]
d={}
d1={}
d2={}
d3={}
with open('essai1.txt') as fic:
	for ligne in fic:
		chaine1 = ligne.split('\t')[:1]
		chaine2= ligne.split('\t')[1:2]
		coco.append(chaine1)
		aster.append(chaine2)
#print aster
 
 
with open('coco_cds.txt') as fichier:
	for line in fichier:
		ch= line.split('\t')[1:2]
		#print ch
		if ch in coco:
			brin=line.split('\t')[3:4]
			if '+'in brin:	
				print "ok"
				cv=line.split('\t')[1:2]
				#print cv
				cv="".join(cv)
				coo=line.strip('\t')
 
				coo=coo.strip('\n')
				liste=coo.split('\t')[7:8]
				seq="".join(liste)
				d[cv]=seq
 
 
for valeur in d.values():
   	valeur=valeur.replace(',','..,..')
	result=valeur.split('..')
	prot=1
	#print result
	nombre1=int(result[0])
	nombre2=int(result[1])
	ok=str((nombre2-nombre1)/3)
	#print ok
	d1[cv]=prot,ok
	#print d1
 
 
num1 = [(cv, ok) for cv, ok in d1.iteritems()]
print num1
 
with open('atro_cds.txt') as fichier:
	for line in fichier:
		ch= line.split('\t')[1:2]
		#print ch
		if ch in aster:
			brin=line.split('\t')[3:4]
			if '+'in brin:	
				print "ok"
				at=line.split('\t')[1:2]
				#print cv
				at="".join(at)
				coo=line.strip('\t')
 
				coo=coo.strip('\n')
				liste=coo.split('\t')[7:8]
				seq="".join(liste)
				d2[at]=seq
				#print d2
 
for valeur in d2.values():
   	valeur=valeur.replace(',','..,..')
	result=valeur.split('..')
	prot2=1
	#print result
	nombre1=int(result[0])
	nombre2=int(result[1])
	ok=str((nombre2-nombre1)/3)
	#print ok
	d3[at]=prot2,ok
	print d3
 
 
num2 = [(at, ok) for cv, ok in d3.iteritems()]
print num2

jai commence ce code mais j'arrive a faire le calcul pour les premiers de chaque liste mais je men mele avec les dico et les listes est ce que vous avez une autre maniere a me dire de faire et en plus il faut que je lui dises que ca fqit 1..4 donc javais commence a mettre prot=1 et donc la jai reussi a mettre dans le dico 1-4 mais il me reste la boucle pour les autres coordonnes et la gestion des virgules qui sont importants

mais je ne sais pas quoi en faire

help please

merci

[wiztricks]

Salut,

Python n'est qu'un langage de programmation.
C'est un bon langage car il permet de mettre rapidement en musique vos idées/algorithmes... Mais qui dit algorithme dit "structurer" un peu ce que vous voulez faire.

Dans votre truc la chose la plus indigeste paraît être la chaîne de caractères:
1412897..1412909,1413103..1413155,1413409..1413531,1413634..1413748,1414228..1414427,1414721..1414891,1415181..1415381,1415548..1415679
Si vous la laissez sous cette forme, le travail que vous allez devoir faire pour la triturer va paralyser vos neurones, pas bon!!!
En fait s'il s'agit d'une suite d'intervalles dont les bornes sont des entiers:
1 - donnez lui un nom: gen_coord
2 - construisez gen_coord à partir de sa représentation "chaine de caractères", i.e. définissez une fonction qui fasse cette transformation.

Code :

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

1
2
3
4
5
6
def get_gen_coord(text):
     l = [] # initialisation
     for e in text.split(','):
           u, v = e.split('..')
           l.append((int(u), int(v)))
     return l

plus pythonique:

Code :

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

get_gen_coord = lambda t: [ (int(u), int(v)) for e in t.split(',') for u,v in [e.split('..')]]

Pour ce qui est de transformer le gen_coord en prot_coord, désolé mais c'est pas très clair, sinon que c'est fonction de la position dans le gen_coord et du quotient de la division par 3:

Code :

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

1
2
3
4
5
def get_prot_coord(gc):
     l = []
     for i, c in enumerate(gc):
         l.append((i+1, (c[1]-c[0])/3)]
     return l

Après, il serait sage de mettre en forme les données lues des différents fichiers. Ce sont des "structures" au sens ou les différents éléments on un nom et éventuellement un type.
En Python (2.7+), vous pouvez assez facilement construire cela avec une chose qui s'appelle namedtuple. Dans votre cas, elle pourrait s'utiliser ainsi:

Code :

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

1
2
from collections import namedtuple
Data = namedtuple('Data', 'scaffold gene_id start dir origine_gene model ident gen_coords')

Les 2 fichiers de CV et ASTER ont la même structure qu'on va pouvoir mettre dans des Data mais il faudra retrouver les entrées correspondantes aux "blasts". Ce qui suggère d'indexer les Data suivant le champ gene_id dans des 'dict' au fur et à mesure qu'on les lit.
=> Nous avons une fonction de lecture des fichier à qui on passe un nom de fichier et que retourne ce qu'il a lu dans un dictionnaire

Code :

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

1
2
3
4
5
6
7
8
9
def readdata(name):
     dd = {}
     with open(name, 'rU') as f:
         for line in f.readlines():
             d = line.split('\t')
             d [-1] = get_gen_coord(d[-1])
             dt = Data(*d)
             dd[dt.gene_id] = dt
     return dd

Maintenant qu'on a un peu simplifié les choses, le coeur du truc devient

Code :

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

1
2
3
4
5
6
7
cv = readdata('fichier de cv')
xx = readdata('fichier de scaffolds')
zz = readdata('ficiher de blast')
for e in zz:
     #on prend le cv et le xx qui vont bien
     #on fait...
     #on affiche le résultat.

Autrement dit dans le corps vous ne mettez que le minimum de verbes qui correspondront aux traitements réalisés dans les différentes fonctions.
- W

[shadow19c]

Bonjour, quand je fais cette fonction

Code :

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

get_gen_coord = lambda t: [ (int(u), int(v)) for e in t.split(',') for u,v in [e.split('..')]]

j'ai une erreure du style

Code :

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

1
2
get_gen_coord = lambda t: [ (int(u), int(v)) for e in t.split(',') for u,v in [e.split('..')]]
ValueError: need more than 1 value to unpack

Mais aussi je ne comprends pas la partie

Code :

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

1
2
3
 
 
xx = readdata('fichier de scaffolds')

car je nai pas de fichier semblable

merci

[wiztricks]

Je n'essaie que de vous donne des idées, c'est vous qui codez.
vos fichiers sont: essai1.txt, coco_cds.txt, astro_cds.txt

get_gen_coord = lambda t: [ (int(u), int(v)) for e in t.split(',') for u,v in [e.split('..')]]

Oubliez cela pour l'instant: il vous sera plus facile de comprendre et de modifier la forme "classique".

le but étant d'écrire quelque chose comme:

Code :

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

1
2
3
4
5
6
coco = readdata ('coco_cds.txt')
astro = readdata ('astro_cds.txt')
essais = readessais('essai1.txt')
for e in essais:
     if e.CV in coco and e.ASTER in astro:
          # on regarde du côté coordonnées.

- W

[shadow19c]

Je vais reprendre ce que jai fais car la je comprends pas tout
en essayant de comprendre petit a petit car la je suis perdu