← Back to branch summary

~ubuntu-branches/debian/sid/python-biopython/sid

« back to all changes in this revision

Viewing changes to Doc/cookbook/motif/motif.tex

  • Committer: Package Import Robot
  • Author(s): Andreas Tille
  • Date: 2014-06-08 10:25:35 UTC
  • mfrom: (1.1.22)
  • Revision ID: package-import@ubuntu.com-20140608102535-ebgnhcpcrum3i3w8
Tags: 1.64+dfsg-1
http://bugs.debian.org/750845
* New upstream version
* Build-Depends raxml only on those architectures where it is available
  Closes: #750845

Show diffs side-by-side

added added

removed removed

Lines of Context:
Doc/cookbook/motif/motif.tex
63
63
This is for now just an empty container, so let's add some sequences to our newly created motif:
64
64
\begin{verbatim}
65
65
from Bio.Seq import Seq
66
 
m.add_instance(Seq("TATAA",m.alphabet))
67
 
m.add_instance(Seq("TATTA",m.alphabet))
68
 
m.add_instance(Seq("TATAA",m.alphabet))
69
 
m.add_instance(Seq("TATAA",m.alphabet))
 
66
m.add_instance(Seq("TATAA", m.alphabet))
 
67
m.add_instance(Seq("TATTA", m.alphabet))
 
68
m.add_instance(Seq("TATAA", m.alphabet))
 
69
m.add_instance(Seq("TATAA", m.alphabet))
70
70
\end{verbatim}
71
71
Now we have a full \verb|Motif| instance, so we can try to get some
72
72
basic information about it. Let's start with length and consensus
139
139
background. We can use the \verb|.log_odds()| method:
140
140
 
141
141
\begin{verbatim}
142
 
 >>> m.log_odds() 
 
142
>>> m.log_odds()
143
143
[{'A': -2.3219280948873622, 
144
144
  'C': -2.3219280948873622, 
145
145
  'T': 1.7655347463629771, 
180
180
(\url{http://jaspar.genereg.net}) stores motifs in both formats, so
181
181
let's look at how we can import JASPAR motifs from instances:
182
182
\begin{verbatim}
183
 
arnt=Motif.read(open("Arnt.sites"),"jaspar-sites")
 
183
arnt = Motif.read(open("Arnt.sites"), "jaspar-sites")
184
184
\end{verbatim}
185
185
and from a count matrix:
186
186
\begin{verbatim}
187
 
srf=Motif.read(open("SRF.pfm"),"jaspar-pfm")
 
187
srf = Motif.read(open("SRF.pfm"), "jaspar-pfm")
188
188
\end{verbatim}
189
189
 
190
190
The \verb|arnt| and \verb|srf| motifs can both do the same things for
230
230
 
231
231
Speaking of exporting, let's look at export functions. We can export to fasta:
232
232
\begin{verbatim}
233
 
>>> print m.format("fasta")
 
233
>>> print(m.format("fasta"))
234
234
> instance 0
235
235
TATAA
236
236
> instance 1
242
242
\end{verbatim}
243
243
or to TRANSFAC-like matrix format (used by some motif processing software)
244
244
\begin{verbatim}
245
 
>>> print m.format("transfac")
 
245
>>> print(m.format("transfac"))
246
246
XX
247
247
TY Motif
248
248
ID 
275
275
The simplest way to find instances, is to look for exact matches of
276
276
the true instances of the motif:
277
277
\begin{verbatim}
278
 
>>> for pos,seq in m.search_instances(test_seq):
279
 
...     print pos,seq.tostring()
 
278
>>> for pos, seq in m.search_instances(test_seq):
 
279
...     print("%i %s" % (pos, seq))
280
280
... 
281
281
10 TATAA
282
282
15 TATAA
284
284
\end{verbatim}
285
285
We can do the same with the reverse complement (to find instances on the complementary strand):
286
286
\begin{verbatim}
287
 
>>> for pos,seq in m.reverse_complement().search_instances(test_seq):
288
 
...     print pos,seq.tostring()
 
287
>>> for pos, seq in m.reverse_complement().search_instances(test_seq):
 
288
...     print("%i %s" % (pos, seq))
289
289
... 
290
290
12 TAATA
291
291
20 TTATA
293
293
 
294
294
It's just as easy to look for positions, giving rise to high log-odds scores against our motif:
295
295
\begin{verbatim}
296
 
>>> for position, score in m.search(test_seq,threshold=5.0):
297
 
...     print position, score
 
296
>>> for pos, score in m.search(test_seq,threshold=5.0):
 
297
...     print("%i %f" % (pos, score))
298
298
... 
299
299
10 8.44065060871
300
300
-12 7.06213898545
345
345
you exactly the same results (for this sequence) as searching for
346
346
instances with balanced threshold with rate of $1000$.
347
347
\begin{verbatim}
348
 
>>> threshold=sd.threshold_balanced(1000)
349
 
>>> for position, score in m.search(test_seq, threshold=threshold):
350
 
...     print position, score
 
348
>>> threshold = sd.threshold_balanced(1000)
 
349
>>> for pos, score in m.search(test_seq, threshold=threshold):
 
350
...     print("%s %s" % (pos, score))
351
351
... 
352
352
10 8.44065060871
353
353
15 8.44065060871
383
383
To show how these functions work, let us first load another motif,
384
384
which is similar to our test motif \verb|m|:
385
385
\begin{verbatim}
386
 
>>> ubx=Motif.read(open("Ubx.pfm"),"jaspar-pfm")
 
386
>>> ubx=Motif.read(open("Ubx.pfm"), "jaspar-pfm")
387
387
<Bio.Motif.Motif.Motif object at 0xc29b90>
388
388
>>> ubx.consensus
389
389
Seq('TAAT', IUPACUnambiguousDNA())
448
448
running the following piece of code:
449
449
 
450
450
\begin{verbatim}
451
 
>>> motifsM = list(Motif.parse(open("meme.out"),"MEME"))
 
451
>>> motifsM = list(Motif.parse(open("meme.out"), "MEME"))
452
452
>>> motifsM
453
453
[<Bio.Motif.MEMEMotif.MEMEMotif object at 0xc356b0>]
454
454
\end{verbatim}
489
489
with the following code:
490
490
 
491
491
\begin{verbatim}
492
 
>>> motifsA=list(Motif.parse(open("alignace.out"),"AlignAce"))
 
492
>>> motifsA = list(Motif.parse(open("alignace.out"), "AlignAce"))
493
493
\end{verbatim}
494
494
 
495
495
Again, your motifs behave as they should:
511
511
shown below (other parameters can be specified as keyword parameters):
512
512
 
513
513
\begin{verbatim}
514
 
>>> command="/opt/bin/AlignACE"
515
 
>>> input_file="test.fa"
516
 
>>> result=Motif.AlignAce(input_file,cmd=command,gcback=0.6,numcols=10)
 
514
>>> command = "/opt/bin/AlignACE"
 
515
>>> input_file = "test.fa"
 
516
>>> result = Motif.AlignAce(input_file, cmd=command, gcback=0.6, numcols=10)
517
517
>>> result
518
518
(<Bio.Application.ApplicationResult instance at 0xf49d78>,
519
519
 <Bio.File.UndoHandle instance at 0xf49d00>,
523
523
Since AlignAce prints all its output to standard output, you can get
524
524
to your motifs by parsing the second member of the result:
525
525
\begin{verbatim}
526
 
motifs=list(Motif.parse(result[1],"AlignAce"))
 
526
motifs = list(Motif.parse(result[1], "AlignAce"))
527
527
\end{verbatim}
528
528
 
529
529