I have a multiple sequence alignment (MSA) of protein sequences in FASTA format. I will like to compute a sequence identity matrix that shows the sequence identity between each sequence and every other sequence in the MSA.
Is there an application/package I can do this with?
A very simple and ultra rapid application works from the Unix/Linux command line with the program Clustal installed. With all your sequences to be aligned in a simple fasta file named "test" you type just one command: clustalw test;clustalw -outputtree=dist -tree -infile=test.aln;l test.aln;l test.dst. This command gives you the following output files:
1, multiple alignment
2, distance matrix
3, a description of a simple phylogenetic (neighbor-joining) tree (test.ph). The latter can be imported or even simply pasted in a program such as FigTree for visualisation
Dear Japheth E. Gado
as Abhijeet Singh alredy wrote, CLustal Omega
https://www.ebi.ac.uk/Tools/msa/clustalo/
can provide you as the identity matrix as output, This matrix could be easly imported in excel of other software.
An example of this is avaialble at the minute 2 of the following video
https://www.youtube.com/watch?v=pujCNbwpUtE
extracted from my blog, ProteoCool
best regards
Manuele
The R package bio3d: http://thegrantlab.org/bio3d/index.html
has the function seqidentity: http://thegrantlab.org/bio3d/reference/seqidentity.html
and approach to plot the PID heatmap with ggplot2 would be:
#calculate PID
mat %
gather(seq2, pid, - 1) %>%
mutate(seq = fct_inorder(seq),
seq2 = fct_inorder(seq2),
pid2 = ifelse(!is.na(pid), paste0(pid * 100, "%"), NA_character_)) %>%
ggplot()+
geom_tile(aes(x = seq, y = seq2, fill = pid)) +
scale_fill_gradient(low = "white",
high = "13ae4b",
na.value = "white",
breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1),
limits = c(0, 1)) +
geom_text(aes(x = seq, y = seq2, label = pid2)) +
theme_classic() +
scale_x_discrete(position = "top") +
theme(legend.position = "none") +
xlab("") +
ylab("") +
coord_equal()
- Badges
- Science topic
- Similar topics
- Filing