Back to index

python-biopython  1.60
simplepredict.py
Go to the documentation of this file.
00001 #!/usr/bin/env python
00002 
00003 # ** What is this code? **
00004 #
00005 # This is an example script for using new biopython modules for
00006 # generating NOE crosspeak peaklists from diagonal peaklist within
00007 # the framework of nmrview.  
00008 # nmrview is not required to run this script, only an installed 
00009 # version of python.
00010 
00011 
00012 # ** What's important? **
00013 #
00014 # The xpktools.py and NOEtools.py modules are what I'm trying to
00015 # demonstrate.  They provide methods and a data class for performing
00016 # some general analysis on NMR data taken directly from peaklist.
00017 # The code in this script will demonstrate how they are used.
00018 
00019 
00020 # ** Who wrote this code? **
00021 # Robert Bussell, Jr.
00022 # rgb2003@med.cornell.edu
00023 
00024 
00025 # ** Running this script **
00026 #
00027 # To run this script on a UNIX/Linux system, make it executable and
00028 # modify the first line of this script to point to python if necessary.
00029 # First try running the code with the peaklist that I provide to get
00030 # the feel of how things work, then you can use your own peaklist if you
00031 # modify the variables under the "INITS" code block to make it work
00032 # with your data.  The modules xpktools and NOEtools can be called from
00033 # your own scripts when you have them in place on your computer.
00034 # NOTE: It is very important to have an intact peaklist.  If you copy and
00035 # paste mine into a file be prepared to remove inappropriate line breaks.
00036 
00037 
00038 # ** Output of this script **
00039 # 
00040 # This script generates a human readable standard output version of the
00041 # NOE coordinates as well as an nmrview peaklist out_example.xpk.
00042 
00043 
00044 
00045 
00046 
00047 
00048 # ***********************************************************************
00049 
00050 # ***** LOAD MODULES *****
00051 import getopt
00052 import string
00053 import sys
00054 
00055 # -- don't need to modify sys.path with the *tools in Biopython
00056 # -- just need Biopython installed somewhere in the PYTHONPATH
00057 #sys.path=[sys.path,"./"]
00058 #sys.path=[sys.path,"/usr/people/robert/development/xpktools"]
00059 from Bio.NMR import xpktools # Contains data classes and functions for .xpk files
00060 from Bio.NMR import NOEtools # A module specific for generate NOE predictions 
00061 
00062 # * * * * * * * * * * MAIN * * * * * * * * * *
00063 
00064 # ***** INITS *****
00065 
00066 inc=1                     # The NOE increment (n where i->i+n and i->i-n are noes)
00067 infn="./noed.xpk"         # Input peaklist
00068 outfn="./out_example.xpk" # Output peaklist
00069 detectatom="H1"           # Directly detected atom
00070 relayatom="N15"           # J-coupling from here to detected atom
00071 fromatom="15N2"           # The other labelled nucleus 
00072 
00073 
00074 
00075 
00076 # *-*-*  First the peaklist is read into a data class from xpktools
00077 # *-*-*  that contains methods for easily extracting information from
00078 # *-*-*  the peaklist file
00079 
00080 peaklist=xpktools.Peaklist(infn) # infn is the name of the xpk file 
00081 
00082 
00083 
00084 
00085 # *-*-* The class attribute residue_dict allows the data lines
00086 # *-*-* to be separated from the header and returned here to
00087 # *-*-* the dictionary <dict> as a list indexed by the assignment
00088 # *-*-* of any of the nuclei in the file -- here, the detected atom
00089 # *-*-* is used 
00090 
00091 dict=peaklist.residue_dict(detectatom)
00092 
00093 
00094 
00095 
00096 # *-*-* As well as the data, the dictionary contains two other entries,
00097 # *-*-* corresponding to the maximum and minimum residues indexed 
00098 
00099 MAXRES=dict["maxres"]
00100 MINRES=dict["minres"]
00101 
00102 
00103 
00104 #****** CALCULATE AND WRITE CROSSPEAK PEAKLIST *****
00105 
00106 
00107 
00108 # *-*-* The peaklist class has a method for writing out the header
00109 # *-*-* information in a format recognizable by nmrview
00110 
00111 peaklist.write_header(outfn) # Write the header to the output file
00112 
00113 
00114 
00115 # *-*-* Predict the i->i+inc and i->i-inc noe positions if possible
00116 # *-*-* Write each one to the output file as they are calculated
00117 
00118 count=0 # A counter that number the output data lines in order
00119 res=MINRES # minimum residue number in the set 
00120 outlist=[] # Holds the output data
00121 
00122 
00123 while (res<=MAXRES):
00124 
00125 
00126 
00127 # *-*-* Predicting the NOE positions based on peak assignment data
00128 # *-*-* is done by supplying the peaklist to and specifying the label
00129 # *-*-* of the origin and detected atom in the NOE transfer as well as
00130 # *-*-* the residues between which the NOE transfer takes places.
00131  
00132   noe1=NOEtools.predictNOE(peaklist,"15N2","H1",res,res+inc)
00133   noe2=NOEtools.predictNOE(peaklist,"15N2","H1",res,res-inc)
00134 
00135 
00136 
00137 # *-*-* The output of predictNOE is in the form of an xpk entry line
00138 # *-*-* suitable for printing to an output file
00139 # *-*-* Additionally, it is possible to extract information easily from
00140 # *-*-* these output lines by using the xpktools.XpkEntry class
00141 
00142   entry1=xpktools.XpkEntry(noe1,peaklist.datalabels)
00143 
00144   if noe1!="":
00145 
00146   # *-*-* Here I'm using the XpkEntry class to gain access to
00147   # *-*-* specific fields in the that make the information
00148   # *-*-* more readable and suitable for creating data tables
00149   # *-*-* This output will be printed to the screen.
00150   # *-*-* The data table contains the assignment, coordinates and
00151   # *-*-* intensity of the resonance.
00152 
00153     print string.split(entry1.fields["15N2.L"],".")[0], "-->", string.split(entry1.fields["N15.L"],".")[0], "\t", entry1.fields["H1.P"], entry1.fields["N15.P"], entry1.fields["15N2.P"],entry1.fields["int"]
00154 
00155 
00156     noe1=noe1+"\012"
00157     noe1=xpktools.replace_entry(noe1,1,count)
00158     outlist.append(noe1)
00159     count += 1
00160 
00161     if noe2!="":
00162       noe2=noe2+"\012"
00163       noe2=xpktools.replace_entry(noe2,1,count)
00164       outlist.append(noe2)
00165       count += 1
00166   res += 1
00167 
00168 # Open the output file and write the data
00169 outfile=open(outfn,'a')
00170 outfile.writelines(outlist) # Write the output lines to the file
00171 outfile.close()