I am doing something extremely basic and stupid. I can call a fasta seq using this
x <- readDNAStringSet(file.choose())
file gets opened, and then I type this
length(x)
[1] 1 But when I use view, I can see nicely that it has 11405 bp. I can use seqinr to get the length, like
length(x[[1]])
But, I can't do anything more, how do I extract the sequence only? And how do I convert it to a DNAString?
•
written
by
c.chakraborty
The readDNAStringSet() reads a set of sequences from a file into an XStringSet object. To see the help page of the XStringSet-class, use:
?XStringSet
There you will find why length(x) returned 1:
length(x): The number of sequences in x.
For an DNAStringSet object, the function you want is width(). If you want a DNAString, you may extract just one element of the DNAStringSet and assign it to a new variable:
xString <- x[[1]]
Now, length(xString) will show the sequence length.
This is great, thanks a bunch! I also found a crude solution. If I load the files by ape and then unlist it, and then ask for length, it works!
Log in to answer this question.
type
width(x). This would print the lengths of all sequences in DNAstringset object.