-4
$\begingroup$

I have a blast file produced. I executed a blast(x) command outputting both "qeseqid" and "sseqid":

QRv313_NP342_d0_h2_l9    YN13213
QRv313_NP9080_d0_h1_l1   YN5345
QRv313_NP123_d0_h1_l7    YN756
QRv313_NP123_d0_h1_l113  YN9768
QRv313_NP654_d0_h2_l6    YN432
QRv313_NP8_d0_h1_l1      YN3242
QRv313_NP756_d0_h1_l2    YN85686

I have written a command in nano within command-line to obtain the following desired output:

NP342    YN13213
NP9080   YN5345
NP123    YN756
NP123    YN9768
NP654    YN432
NP8_d0   YN3242
NP756    YN85686

I have written a nano script to provide me a tab delimited column of my query and subject id. I am just having trouble moving forward from here. I am unsure as to how I would modify my script to provide me with my desired output.

import sys
file_object = open(sys.argv[1])

for my_data in file_object:

  list =  my_data.split("\t")

  print (list [0], list [1])

Is there a way t alter my command so I can receive the desired output?

Any suggestions would be kindly appreciated!

$\endgroup$
3
$\begingroup$

as far as I understand your question, you want to reformat your output file that consists of two columns:

QRv313_NP342_d0_h2_l9    YN13213

to extract the second "field" from your first column, where fields are delimited by _:

NP342    YN13213

In general, your code is close. You could just split again by _:

import sys
file_object = open(sys.argv[1])

for my_data in file_object:
  list =  my_data.strip().split("\t")
  print (list[0].split("_")[1], list[1], sep="\t")

Now when I teach students, I usually tell them to always expect the unexpected, so you might also want to add some minimal error checking, so here's a version that checks some minimal assumptions:


from sys import argv
from os import path

usage = f"""
usage: {argv[0]} <blast_output>
"""


def main(in_file):
    with open(argv[1], "r") as in_handle:
        for lno, line in enumerate(in_handle):
            line = line.strip()
            cols = line.split("\t")
            if len(cols) != 2:
                raise ValueError(
                    f"Error in line {lno}: Expected two columns,"
                    " found {len(cols)}. Offending line: {cols}"
                )
            query = cols[0].split("_")
            if len(query) != 5:
                raise ValueError(
                    f"Error in line {lno}: Expected five fields"
                    "in query string, found {len(query)}. Offending line {query}"
                )
            print(query[1], cols[1], sep="\t")


if __name__ == "__main__":
    if len(argv) != 2:
        print(usage)
        raise ValueError("Illegal number of arguments")
    if not path.exists(argv[1]):
        raise FileNotFoundError(f"File {argv[1]} could not be found. Check spelling")

    main(argv[1])

And this is what it produces:

>python3 ex.py in.txt
NP342   YN13213
NP9080  YN5345
NP123   YN756
NP123   YN9768
NP654   YN432
NP8     YN3242
NP756   YN85686
$\endgroup$
2
  • 1
    $\begingroup$ Yup, that was a typo. The keyword is meant to be sep="\t". Corrected in the script $\endgroup$ Nov 29 '20 at 20:32
  • $\begingroup$ Keep track of the ID's you've encountered within python and perform any filtering you want according to that. Basic python programming is outside the scope of this site. $\endgroup$
    – Devon Ryan
    Nov 30 '20 at 11:25

Not the answer you're looking for? Browse other questions tagged or ask your own question.