Forum

Nome Utente:
Password:
Riconoscimi automaticamente
 Tutti i Forum
 MolecularLab
 Bioinformatica
 ricerca siti di restrizione con biopython
 Nuova Discussione  Nuovo Sondaggio Nuovo Sondaggio
 Rispondi Aggiungi ai Preferiti Aggiungi ai Preferiti
Cerca nelle discussioni
I seguenti utenti stanno leggendo questo Forum Qui c'č:

Aggiungi Tag Aggiungi i tag

Quanto č utile/interessante questa discussione:

Autore Discussione  

HnACP
Nuovo Arrivato


Prov.: MI
Cittā: Milano


78 Messaggi

Inserito il - 13 maggio 2009 : 09:46:07  Mostra Profilo  Visita l'Homepage di HnACP Invia a HnACP un Messaggio Privato  Rispondi Quotando
ciao a tutti!
ho scritto un codice per trovare siti di restrizione in una sequenza fasta, ma nn riesco a dare l'output, nel senso che quando lo faccio girare nn mi da niente. Copio il codice qua: qualcuno mi sa aiutare?
Se funziona magari puo' essere di aiuto a qualcuno!

from Bio import SeqIO
from string import *
handle = open("unr.fasta")
for sequnr in SeqIO.parse(handle, "fasta") :
 
 def digest(sequnr, enz) :
  EcoRI = 'GAATTC'
  BamHI = 'GGATCC'
  HindIII = 'AAGCTT'
  enz = [ EcoRI, BamHI, HindIII ]
    
  res = []

  site = sequnr.find(enz)
  while site != -1 :
       res.append(site)
       site = sequnr.find(enz, site + 1)
  return res

handle.close()

dallolio_gm
Moderatore


Prov.: Bo!
Cittā: Barcelona/Bologna


2445 Messaggi

Inserito il - 13 maggio 2009 : 10:04:51  Mostra Profilo  Visita l'Homepage di dallolio_gm  Clicca per vedere l'indirizzo MSN di dallolio_gm Invia a dallolio_gm un Messaggio Privato  Rispondi Quotando
ciao,
forse l'indentazione é sbagliata, o si é persa qualche riga nel post, perķ non capisco bene se stai definendo la funzione digest per ogni esecuzione del ciclo for.
Se il codice é quello, non ti funziona perché non stai chiamando la funzione digest(), la stai solo definendo.

Innanzitutto, una buona regola quando si scrivono programmi é quella di pensare subito a un test case per verificare che i risultati siano corretti. Per esempio, dovresti scrivere un paio di sequenze casuali e verificare che il tuo script funzioni correttamente.

Secondo, su biopython esiste un modulo chiamato Restriction, per lavorare su enzimi di restrizione. Perķ il modulo sembra essere piuttosto vecchio e almeno a me, non funziona.
Ti sconsiglio di utilizzarlo, e sarebbe il caso di avvertire gli sviluppatori (posso farlo io se vuoi).

Terzo, per lavorare su motivi e matrici pssm, ti consiglio questo modulo:
- http://fraenkel.mit.edu/TAMO/TAMO_tutorial.html
E' ottimo e se dai una occhiata al tutorial, puoi vedere che esiste un metodo .scan() che permette di cercare le occorrenze del motivo su una sequenza.
In ogni caso, ti sconsiglio di utilizzare il semplice .find, piuttosto dai una occhiata alle espressioni regolari (modulo re).

ok, adesso aspetto una tua risposta :)

Il mio blog di bioinformatics (inglese): BioinfoBlog
Sono un po' lento a rispondere, posso tardare anche qualche giorno... ma abbiate fede! :-)
Torna all'inizio della Pagina

dallolio_gm
Moderatore


Prov.: Bo!
Cittā: Barcelona/Bologna


2445 Messaggi

Inserito il - 13 maggio 2009 : 10:24:13  Mostra Profilo  Visita l'Homepage di dallolio_gm  Clicca per vedere l'indirizzo MSN di dallolio_gm Invia a dallolio_gm un Messaggio Privato  Rispondi Quotando

from Bio import SeqIO
from TAMO import MotifTools

# from string import * # non necessario

handle = open("unr.fasta", 'r')

EcoRI = MotifTools.Motif_from_text('GAATTC')
BamHI = MotifTools.Motif_from_text('GGATCC')
HindIII = MotifTools.Motif_from_text('AAGCTT')

enzimes = [ EcoRI, BamHI, HindIII ]

for sequnr in SeqIO.parse(handle, "fasta"):
  for enz in enzimes:
    digest(sequnr, enz)


def digest(sequnr, enz) :
  """
  trova tutte le occorrenze di un motivo in una sequenza. (corretto?)

  esempio:
  >>> from StringIO import StringIO
  >>> from TAMO import MotifTools
  >>> from Bio import SeqIO
  >>> fastafile = StringIO('>seq1\nGGATcCCGGATGCATGCTA')
  >>> fastafile.seek(0)
  >>> enz = MotifTools.Motif_from_text('GGATCC')

  >>> print digest(SeqIO.parse(fastafile, 'fasta'), enz)
  (['GGATcC'], [(0, 5)], [11.685195480635185])
  """    

  occurrencies = enz.scan(sequnr.seq.tostring())
  return occurencies

Probabilmente puoi spostare alcune chiamate all'interno della funzione, come la creazione del motivo e il parsing del file.
Scrivi prima il doctest e poi ti sará facile decidere come impostarlo.

ultima cosa: non usare il metodo find di biopython, non é case insensitive e non distingue le N come basi speciali.

Il mio blog di bioinformatics (inglese): BioinfoBlog
Sono un po' lento a rispondere, posso tardare anche qualche giorno... ma abbiate fede! :-)
Torna all'inizio della Pagina

HnACP
Nuovo Arrivato


Prov.: MI
Cittā: Milano


78 Messaggi

Inserito il - 13 maggio 2009 : 23:12:26  Mostra Profilo  Visita l'Homepage di HnACP Invia a HnACP un Messaggio Privato  Rispondi Quotando
ciao!
quel codice, un po arrangiato da me, l'avevo trovato sul tutorial di biopython dell'istituto pasteaur, solo che č cosi stringato che nn si capisce niente...

in pratica volevo, da una sequenza in formato fasta, vedere se dei dati siti di restrizione erano presenti in quella sequenza e quanti erano. Con un altro script sono riuscito solo che mi indicava solo il primo sito e poi si fermava! o cmq se nn sono siti di restrizione magari cercare una data sequenza di un promotore, ma penso che la procedura sia la stessa.

quel modulo TAMO sembra molto interesseante, solo che su linux nn riesco ad installarlo...

sul tutorial usa quasi sempre find... che poi ha funzionato solo una volta con me..

alla fine io mi basavo solo sul tutorial ufficiale di biopython e quello francese.. ma vedo che hanno cose datate... potrebbero aggiornarsi un po, no?

visto che nn riesco ad installare il modulo TAMO nn c'č il modo di fare lo script in puro python?
Torna all'inizio della Pagina

dallolio_gm
Moderatore


Prov.: Bo!
Cittā: Barcelona/Bologna


2445 Messaggi

Inserito il - 14 maggio 2009 : 00:36:02  Mostra Profilo  Visita l'Homepage di dallolio_gm  Clicca per vedere l'indirizzo MSN di dallolio_gm Invia a dallolio_gm un Messaggio Privato  Rispondi Quotando
Che problemi hai ad installare TAMO?
Basta scompattare il file e dare 'sudo python setup.py install'.
E' richiesto il modulo numpy... lo puoi installare tramite easy_install, o come pacchetto nella tua distro.

Se non vuoi utilizzare moduli esterni (male... mai riscrivere quello che giā esiste), usa le espressioni regolari.
Qualcosa del genere:

>>> import re
>>> seq1 = 'ACGATGGCTACGTACGTACGTACGA'
>>> seq2 = 'AAAAAAAACCAAAAAAAAAAAAAAA'

>>> re.findall('AC', seq1)
['AC']

>>> for match in re.finditer('AA', seq2):
        print match.start(), match.end(), match.group()

# nota come anche con le regex il comportamento é incorretto, perché non si tiene conto che una N rappresenta qualsiasi nucleotide:
>>> seq3 = 'ANAA'
>>> re.findall('AA', seq3)

Il mio blog di bioinformatics (inglese): BioinfoBlog
Sono un po' lento a rispondere, posso tardare anche qualche giorno... ma abbiate fede! :-)
Torna all'inizio della Pagina

HnACP
Nuovo Arrivato


Prov.: MI
Cittā: Milano


78 Messaggi

Inserito il - 14 maggio 2009 : 23:05:51  Mostra Profilo  Visita l'Homepage di HnACP Invia a HnACP un Messaggio Privato  Rispondi Quotando
ciao! si ho provato ad installarlo ma ad un certo punto mi da errore:

creating build/temp.linux-i686-2.6/TAMO/MD/MDsupport_source
gcc -pthread -fno-strict-aliasing -DNDEBUG -g -fwrapv -O2 -Wall -Wstrict-prototypes -fPIC -I/usr/include/python2.6 -c TAMO/MD/MDsupport_source/MDsupport.cxx -o build/temp.linux-i686-2.6/TAMO/MD/MDsupport_source/MDsupport.o
gcc: error trying to exec 'cc1plus': execvp: No such file or directory
error: command 'gcc' failed with exit status 1
Torna all'inizio della Pagina
  Discussione  

Quanto č utile/interessante questa discussione:

 Nuova Discussione  Nuovo Sondaggio Nuovo Sondaggio
 Rispondi Aggiungi ai Preferiti Aggiungi ai Preferiti
Cerca nelle discussioni
Vai a:
MolecularLab.it © 2003-18 MolecularLab.it Torna all'inizio della Pagina