Remove spaces from fasta file in python

Question

[anasjamshed]

I have downloaded Severe acute respiratory syndrome coronavirus 2 isolate FASTA sequence but when I tried to get the exact length and C and G count I got the wrong answer every time.

Script:

genome = open("C://Users//USER//Desktop/Corona.fasta")
dna = genome.read()
dna1 = dna.rstrip("\n")
print(len(dna1))

The exact length of the sequence is 29903 without a header 30018 but I don't understand why I am not getting 30018 although I am using strip function as well

Plz help me

[cpad0112]

Please use biopython or other fasta parsing libraries, for reading fasta files and subsequent operations. Post example file for better understanding the query.

[anasjamshed]

I have already solved it by using biopython but I want to do it without using any library?

[cpad0112]

seq=[]
with open ("test.fa","r") as w:
    for i in w:
        if i.startswith (">"):
            _head=i.strip("\n>")
            _seq= next(w, '').strip('').replace(" ","")
            gc_count=sum(map(_seq.upper().count,["G","C"]))
            gc_perc=gc_count/len(_seq)
            _temp_seq=[_head, len(_seq),gc_count,gc_perc]
            seq.append(_temp_seq)

[print (i) for i in seq] 

Check this code. This prints fasta header, length of the sequence, Number of G and C, % GC. It strips all the spaces in the seqence. This code expects the sequence in a single line in a fasta record

[anasjamshed]

what is with?

[anasjamshed]

['gi|1798174254|ref|NC_045512.2| Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, complete genome', 71, 25, 35.2112676056338]

It is also producing wrong length and GC content

[cpad0112]

You have missed this part of my post: This code expects the sequence in a single line in a fasta record. From your post, I can see that fasta record used is multiline record.

Following is the comparison of results from fairly known tool and the script and I have used this sequence for test:

$ python test.py 
['NC_045512.2 Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, complete genome', 29904, 11355, 0.379715088282504]

$ seqkit fx2tab -n -g seq.fa 
NC_045512.2 Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, complete genome 37.97

GC content from the script : 0.379715088282504 and GC content from the tool: 37.97. Code doesn't multiply GC content by 100.

[anasjamshed]

which tool you are using?

[jared.andrews07]

It says quite clearly in the command - seqkit.

Answer 1

[jared.andrews07]

There are numerous examples of substring counting that can be found with a bit of googling, e.g. here.

Please post a link to the sequence/file in question if you want specific help with your issue. My guess is that it could be a platform-specific issue (try dna.rstrip("\n\r"), which should remove all newline characters. Alternatively, just plan dna.rstrip() should also remove all newline and whitespace characters (spaces, tabs).

You are reading the whole file into a single string, from which only the last newline character is being stripped. It is easier to replace all instances of the newline character since it is going to be interspersed throughout your read.

dna1 = dna.replace('\n', '')

[anasjamshed]

rstrip just removing 2 spaces

[jared.andrews07]

Oh, you are reading in the entire file at once, so my answer will indeed only remove the last newline character. I will edit my answer.

[anasjamshed]

how I remove it?

[anasjamshed]

how can I remove the header?

[jared.andrews07]

You can just read the file line by line while trimming and skip the header:

dna = []
with open("C://Users//USER//Desktop/Corona.fasta") as genome:
    # Skip header.
    genome.readline()

    # Iterate through each line, strip newline, and add to your dna list.
    for line in genome:
        line = line.strip()
        dna.append(line)

# Concatenate all elements of the list to a string with no space between them.
dna1 = "".join(dna)
print(len(dna1))