tankdthedruid

NCBI_record_parser.py

Jan 22nd, 2012
103
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Python 7.24 KB | None | 0 0
  1. #   ~/NCBI_record_parser.py
  2. #
  3. #FILE_INPUT = open("blast_output.xml")
  4. #MINIMUM_INT = 20
  5. #from Bio.Blast import NCBIXML
  6. #PARSE_FILE = NCBIXML.parse(FILE_INPUT)
  7. #READ_LINE = PARSE_FILE.next()
  8. #
  9. #
  10. #
  11. #for READ_LINE in PARSE_FILE:
  12. #   for alignment in READ_LINE.alignments:
  13. #       for hsp in alignment.hsps:
  14. #           if hsp.identities > MINIMUM_INT:
  15. #               print \n***Seperator*** \n\n ALIGNMENT_TITLE: , alignment.title
  16. #               print length: , alignment.length
  17. #               print HSP_SCORE: , hsp.score
  18. #               print HSP_BITS: , hsp.bits
  19. #               print HSP_EXPECT: , hsp.expect
  20. #               print HSP_NUM_ALIGNMENTS: , hsp.num_alignments
  21. #               print HSP_IDENTITIES: , hsp.identities
  22. #               print HSP_POSITIVES: , hsp.positives
  23. #               print HSP_GAPS: , hsp.gaps
  24. #               print HSP_STRAND: , hsp.strand
  25. #               print HSP_FRAME: , hsp.frame
  26. #               print HSP_QUERY: , hsp.query[0:75] + ...
  27. #               print HSP_: QUERY_START, hsp.query_start
  28. #               print HSP_MATCH: , hsp.match[0:75] + ...
  29. #               print HSP_SBJCT: , hsp.sbjct[0:75] + ...
  30. #               print HSP_SBJCT_START: , hsp.sbjct_start
  31. #in the format below i wrote the variable names in all upper case thought it might be a little easier to distinguish at a glance.
  32. #
  33. #The following variables are declared below:
  34. #
  35. #      --- FILE_INPUT
  36. #      --- MINIMUM_INT
  37. #      --- VARIABLE_ID
  38. #      --- PARSE_FILE
  39. #      --- READ_LINE
  40. #  
  41. #
  42. #
  43. #   set the variable FILE_INPUT equal to the data contained in the file blast_output.xml.
  44. FILE_INPUT = open("blast_output.xml")
  45. #
  46. #   established the minimum threshold integer to 20.
  47. MINIMUM_INT = 20
  48. #
  49. #   from BioPython import the module NCBIXML
  50. from Bio.Blast import NCBIXML
  51. #
  52. #   in the string below the parse() function contained in the module NCBIXML was called to run against the data contained in the variable FILE_INPUT.
  53. PARSE_FILE = NCBIXML.parse(FILE_INPUT)
  54. #
  55. #   now we take the string we just created above and we are going to add the next() function. Thie next function tells the interpreter that when parsing the input, move just one line at a time.
  56. READ_LINE = PARSE_FILE.next()
  57. #
  58. #
  59. #   so a horrible visual for this might look like:
  60. #
  61. #
  62. #      --- NCBIXML (Module)  >  parse()  >  next()
  63. #
  64. #
  65. #   its like we are gradually narrowing down the amount of data that will be look at by the interpreter at any given time.
  66. #
  67. #
  68. #      --- parse() == read entire file at a time
  69. #  
  70. #      --- next() == read a single line from a file at a time
  71. #
  72. #
  73. #
  74. #   now we begin a for loop to go through our data line by line and process the data. So since above we essentially made READ_LINE a parsing function, when it gets called below, it is being told to run against the data contained in variable PARSE_FILE.
  75. for READ_LINE in PARSE_FILE:
  76. #
  77. #
  78. #   so since:
  79. #
  80. #      --- FILE_INPUT = open("blast_output.xml")
  81. #  
  82. #          --- A.K.A FILE_INPUT is equal to the data contained in blast_output.xml
  83. #
  84. #      --- PARSE_FILE = NCBIXML.parse(FILE_INPUT)
  85. #  
  86. #          --- A.K.A PARSE_FILE is equal to the data contained in blast_output.xml and the parse() function is like it stating it is going to stream edit the data in some way.
  87. #  
  88. #      --- READ_LINE = PARSE_FILE.next()
  89. #  
  90. #          --- A.K.A READ_LINE is equal to the data from blast_output.xml, the intention to stream edit, and finally the next() function was chosen as the parsing method. next() is a line-by-line parser
  91. #
  92. #   now this is where things can get confusing.
  93. #
  94. #   the command begins what is called commonly refered to in programming as nesting. this would be considered a for loop nest. what happens in for loop nesting can sometimes be difficult to follow.
  95. #
  96. #   i am going try to break this loop down, by explaining each piece. In the command above:
  97. #
  98. #
  99. #      --- for READ_LINE in PARSE_FILE:     ;       //  1st loop
  100. #  
  101. #   we know its parsing line by line. however when the line below is nested inside:
  102. #
  103. #      --- for alignment in READ_LINE.alignments:       ;       //  2nd loop
  104. #  
  105. #   it says to the interpreter, "parse the whole file, one line at a time, but at the end of each line, i want you to pause and run the command below."
  106. #
  107. #   so then the second loop begins. this one was a little tricky for me too figure out until i did some research. While reading the biopython documentation I ran across this pdf document of the biopython cookbook.
  108. #
  109. #      --- http://biopython.org/DIST/docs/tutorial/Tutorial.pdf     ;       //  Biopython Tutorial and Cookbook Link
  110. #  
  111. #   there is a section that begins on page 87 of that document called Parsing BLAST output. This chapter really helped put things into focus for me. the most important aspect of this chapter to me was the illustrated database on pages 91-92.
  112. #
  113. #   what had me the most confused was that i did not realize that the alignment variables used below are set globally when we import the NCBIXML module. additonionally, even though the variables look almost identical, the data they contain references two different tiers of there record database.
  114. #
  115. #      --- alignment | alignments       ;       //      only one letter difference try not to confuse the two.
  116. #
  117. #   so the below loop could be read something like this:
  118. #
  119. #      --- read the document blast_output.xml, at the end of each line pause, then start to parse the alignment entries in the database. this part of the loop parses out individual records at a time.
  120. #  
  121. #   remember nesting is all about narrowing your search results
  122. #
  123. #
  124.     for alignment in READ_LINE.alignments:
  125. #
  126. #
  127. #   quick recap of wher this program has brought us:
  128. #
  129. #FILE_INPUT = open("blast_output.xml")      ;       //  Open file blast_output.xml
  130. #MINIMUM_INT = 20       ;       //  Set minimum Threshold value
  131. #from Bio.Blast import NCBIXML      ;       //  import the biopython module NCBIXML
  132. #PARSE_FILE = NCBIXML.parse(FILE_INPUT)     ;       //  define PARSE_FILE as the parsed output of blast_output.xml
  133. #READ_LINE = PARSE_FILE.next()      ;       //  define line-by-line parsing of the xml file by using the builtiin parser tools next() function. read file one line at a time.
  134. #for alignment in READ_LINE.alignments:     ;       //  at the end of each line, loop through and parse the alignment entries into their own individual records.
  135. #for hsp in alignment.hsps:     ;       //  each time an individual record is read above, take that record and parse it. this final for loop parses the records into their smallest form, individual strings of data.
  136. #if hsp.identities > MINIMUM_INT:       ;       //  finally parse each string holds an identities value, and compare that number to the threshold. if the value of the identities string is greater than the threshold integer, continue executing this program. otherwise ignore the string and continue the loop until completion.
  137. #
  138. #
  139.         for hsp in alignment.hsps:
  140.             if hsp.identities > MINIMUM_INT:
  141.                 print '\n***Seperator*** \n\n ALIGNMENT_TITLE: ', alignment.title
  142.                 print 'length: ', alignment.length
  143.                 print 'HSP_SCORE: ', hsp.score
  144.                 print 'HSP_BITS: ', hsp.bits
  145.                 print 'HSP_EXPECT: ', hsp.expect
  146.                 print 'HSP_NUM_ALIGNMENTS: ', hsp.num_alignments
  147.                 print 'HSP_IDENTITIES: ', hsp.identities
  148.                 print 'HSP_POSITIVES: ', hsp.positives
  149.                 print 'HSP_GAPS: ', hsp.gaps
  150.                 print 'HSP_STRAND: ', hsp.strand
  151.                 print 'HSP_FRAME: ', hsp.frame
  152.                 print 'HSP_QUERY: ', hsp.query[0:75] + '...'
  153.                 print 'HSP_: QUERY_START', hsp.query_start
  154.                 print 'HSP_MATCH: ', hsp.match[0:75] + '...'
  155.                 print 'HSP_SBJCT: ', hsp.sbjct[0:75] + '...'
  156.                 print 'HSP_SBJCT_START: ', hsp.sbjct_start
Advertisement
Add Comment
Please, Sign In to add comment