3

I am trying to compare two gene sequences:

sequence_1 <- "MPHLENVVLCRESQVSILQSLFGERHHFSFPSIFIYGHTASGKTYVTQTLLKTLELPHVFVNCVECFTLRLLLEQILNKLNHLSSSEDGCSTEITCETFNDFVRLFKQVTTAENLKDQTVYIVLDKAEYLRDMEANLLPGFLRLQELADRNVTVLFLSEIVWEKFRPNTGCFEPFVLYFPDYSIGNLQKILSHDHPPEYSADFYAAYINILLGVFYTVCRDLKELRHLAVLNFPKYCEPVVKGEASERDTRKLWRNIEPHLKKAMQTVYLREISSSQWEKLQKDDTDPGQLKGLSAHTHVELPYYSKFILIAAYLASYNPARTDKRFFLKHHGKIKKTNFLKKHEKTSNHLLGPKPFPLDRLLAILYSIVDSRVAPTANIFSQITSLVTLQLLTLVGHDDQLDGPKYKCTVSLDFIRAIARTVNFDIIKYLYDFL"

sequence_2 <- "MEEEAPRFNVLEEAFNGNGNGCANVEATQSAILKVLTRVNRFQMRVRKHIEDNYTEFLPNNTSPDIFLEESGSLNREIHDMLENLGSEGLDALDEANVKMAGNGRQLREILLGLGVSEHVLRIDELFQCVEEAKATKDYLVLLDLVGRLRAFIYGDDSVDGDAQVATPEVRRIFKALECYETIKVKYHVQAYMLQQSLQERFDRLVQLQCKSFPTSRCVTLQVSRDQTQLQDIVQALFQEPYNPARLCEFLLDNCIEPVIMRPVMADYSEEADGGTYVRLSLSYATKEPSSAHVRPNYKQVLENLRLLLHTLAGINCSVSRDQHVFGIIGDHVKDKMLKLLVDECLIPAVPESTEEYQTSTLCEDVAQLEQLLVDSFIINPEQDRALGQFVEKYETYYRNRMYRRVLETAREIIQRDLQDMVLVAPNNHSAEVANDPFLFPRCMISKSAQDFVKLMDRILRQPTDKLGDQEADPIAGVISIMLHTYINEVPKVHRKLLESIPQQAVLFHNNCMFFTHWVAQHANKGIESLAALAKTLQATGQQHFRVQVDYQSSILMGIMQEFEFESTHTLGSGPLKLVRQCLRQLELLKNVWANVLPETVYNATFCELINTFVAELIRRVFTLRDISAQMACELSDLIDVVLQRAPTLFREPNEVVQVLSWLKLQQLKAMLNASLMEITELWGDGVGPLTASYKSDEIKHLIRALFQDTDWRAKAITQIV"

using the align_local function from the textreuse package. My input is:

library(textreuse)
align_local(sequence_1, sequence_2)

and I get the error:

Error in b_out[out_i] <- b_orig[row_i - 1] : replacement has length zero
In addition: Warning message:
Multiple optimal local alignments found; selecting only one of them.

I've tried tinkering with the alignment score and the mismatch score, but to no avail. Any advice would be appreciated.

MrFlick
  • 195,160
  • 17
  • 277
  • 295
Murph
  • 33
  • 2

2 Answers2

4

The textreuse package is intended for natural language. Under no circumstances should you use it for aligning gene sequences. (I am the package author.) You probably want the Biostrings package from Bioconductor.

The problem is that the align_local() function expects there to be multiple words, as indicated by spaces or punctuation, because it aligns word by word not character by character. The function would work if you put spaces between the bases in your gene sequence. But I'm not going to explain how to do that because, again, you should not be using a natural language package for aligning genes.

Lincoln Mullen
  • 6,257
  • 4
  • 27
  • 30
  • I have a question @Mullen: How about if I want to use a character tokenizer in textreuse? is this possible? I'd like to form my shingles on a character level. – nhern121 Feb 21 '18 at 19:57
  • @Nestorghh When you create a corpus using `TextReuseCorpus` you can pass a function to the `tokenizers` argument to tokenize the texts however you wish. – Lincoln Mullen Feb 26 '18 at 14:45
  • @Mullen, thanks! I just created a tokenizer based on character but I'm getting this warning : ``Skipping document with ID 'r908' because it has too few words to create at least two n-grams with n = 3`` . but ``r908`` document has 10 characters. I'm dealing with short text : minimum length is 3, so with ``k=3`` (3 character shingle) should work but keep getting this. However, if I apply my function out of ``TextReuseCorpus``, it works and return the expected output. Any suggestion? – nhern121 Feb 26 '18 at 18:54
1

The problem here is that lsh_compare function from textreuse package is meant for analyzing text documents and detecting passages which have been re-used. Which means that it works with spaced out words in a sentence.

My suggestion would be to try to find a package that is more suitable for handling genes.

e.g. dotPlot function from seqinr gives you a visual representation of the comparison.

Marija T.
  • 53
  • 1
  • 7