Get sequences command give me an error

[rpestana94]

I'm trying to read a few fasta files which various entries, a few of them have symbols like J or X or N and this is giving me an error, how I can read this kind of files.

The code I'm using is:

import biotite.sequence as seq
import biotite.sequence.io.fasta as fasta

fasta_file = fasta.FastaFile.read("test_few_sample.fasta")
all_sequences = fasta.get_sequences(fasta_file)

And the error I'm getting is:
Traceback (most recent call last):
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/io/fasta/convert.py", line 218, in _convert_to_sequence
prot_seq = ProteinSequence(seq_str)
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/seqtypes.py", line 472, in init
super().init(sequence)
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/sequence.py", line 147, in init
self.symbols = sequence
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/sequence.py", line 183, in symbols
self._seq_code = alph.encode_multiple(value, dtype)
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/alphabet.py", line 380, in encode_multiple
return encode_chars(alphabet=self._symbols, symbols=symbols)
File "src/biotite/sequence/codec.pyx", line 74, in biotite.sequence.codec.encode_chars
biotite.sequence.AlphabetError: Symbol 'J' is not in the alphabet

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
File "matrix.py", line 18, in
all_sequences = fasta.get_sequences(fasta_file)
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/io/fasta/convert.py", line 93, in get_sequences
seq_dict[header] = _convert_to_sequence(seq_str)
File "/home/arcadio/anaconda3/envs/py38/lib/python3.8/site-packages/biotite/sequence/io/fasta/convert.py", line 228, in _convert_to_sequence
raise ValueError("FASTA data cannot be converted either to "
ValueError: FASTA data cannot be converted either to 'NucleotideSequence' nor to 'ProteinSequence'

Thanks in advance.

[padix-key]

'J' is currently not a symbol in the amino acid alphabet. Hence neither a NucleotideSequence or ProteinSequence can be created from the sequences in your FASTA file. There are two possible solutions to this issue, both using the low-level API of FastaFile that returns strings instead of Sequence objects:

1. Read the sequences as string, replace the symbol with an appropriate replacement and create a ProteinSequence

sequences = {header : ProteinSequence(seq_str.replace("J", "L")) for header, seq_str in fasta_file.items()}

2. Read the sequences as string and create GeneralSequence objects with a custom alphabet

alphabet = LetterAlphabet(ProteinSequence.alphabet.get_symbols() + ["J"])
sequences = {header : GeneralSequence(alphabet, seq_str) for header, seq_str in fasta_file.items()}

[rpestana94]

Is there a way to convert a GeneralSequence to a ProteinSequence object? I'm trying to do a multiple alignments with clustal and give me an error saying that "The software cannot align sequences of type GeneralSequence: No support for custom substitution matrices" but running by console and don't have any problem, is there is a way to fix that without changing the J for a L?

[padix-key]

I fear this is not possible in the current implementation.

Is replacing an issue in this case? As far as I know Clustal Omega uses the classic Gonnet or BLOSUM substitution matrices, which do not support 'J' as symbol. So I assume that Clustal Omega internally converts 'J' to some other symbol anyway for the alignment.

[rpestana94]

Is there a way to read an alignment file as an alignment object?

[padix-key]

There is biotite.sequence.io.fasta.get_alignment() (https://www.biotite-python.org/apidoc/biotite.sequence.io.fasta.get_alignment.html), but this won't work for sequences with 'J' for the same reasons as above.
However, you can use Alignment.trace_from_strings(seq_strings) to build an Alignment yourself:

alphabet = LetterAlphabet(ProteinSequence.alphabet.get_symbols() + ["J"])
seq_strings = list(fasta_file.values())
# Remove gaps for creation of sequences
sequences = [GeneralSequence(alphabet, seq_str.replace("-", ""))
             for seq_str in seq_strings]
trace = Alignment.trace_from_strings(seq_strings)
return Alignment(sequences, trace, score=None)

I have not tested this snippet, but I hope this code will work nonetheless.

[rpestana94]

Thanks, I actually deleted the sequences with a J to save the trouble, there were just 13 from 156000 approx so that not affect too much