This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Extracting N positions from fasta file

Hi I am very new to this so appologies if this is a simple question - I have been trying to figure this out for days to no avail and my python skills are not quite there yet!

I am trying to extract all N positions from novel bacterial sequences which have been aligned to a member of the same genus. I would like the start and end positions of all N motifs eg. GTCAGNNNNNTGGT

Is there an existing tool / how could I go about creating this in python?

Many thanks.

sequence

Do you want to search for a specific pattern or for every location in which an 'N' is present?

every location at which a N is present. I have something like 36 different strains and would like to produce a list of N locations in each FASTA file and then compare these lists to find the unique N locations for each strain.

So the output would be the chromosomal locations, right? Sounds like a job that can be done using Biopython. What have you tried?

If the sequence is not long you can do it without software, open the file fasta by the wordpad

Just use SeqKit. shenwei356 has even provided a detailed example below.

2 answers

Use SeqKit. SeqKit supports Windows/Linux/Mac OS X.

$ echo -en '>seq\nGTCAGNNNNNTGGT\n' | seqkit locate --ignore-case --only-positive-strand --pattern "N+" 
seqID   patternName     pattern strand  start   end     matched
seq     N+      N+      +       6       10      NNNNN

aligned result

$ echo -en '>seq\nGTCAGNNNNNTGGT\n' | seqkit locate --ignore-case --only-positive-strand --pattern "N+" | column -t
seqID  patternName  pattern  strand  start  end  matched
seq    N+           N+       +       6      10   NNNNN

or from file and save to file

$ seqkit locate --ignore-case --only-positive-strand --pattern "N+" seqs.fa > result.xls

picard ScatterIntervalsByNs https://broadinstitute.github.io/picard/command-line-overview.html#ScatterIntervalsByNs

Writes an interval list based on splitting a reference by Ns. This tool identifies positions in a reference where the bases are 'no-calls' and writes out an interval-list using the resulting coordinates. This can be used to create an interval list for whole genome sequence (WGS) for e.g. scatter-gather purposes, as an alternative to using fixed-length intervals. The number of contiguous nocalls that can be tolerated before creating a break is adjustable from the command line.

Log in to answer this question.