I have run OrthoFinder on a set of 11 genomes and got the results in a folder. From the output folder, I saved all the single-copy gene sequences in a single FASTA file and now wish to remove all gaps using Gblocks before using them for constructing a phylogenomic tree. But I don't know the command line used for running Gblocks. Any help will be appreciated!
4
-
$\begingroup$ I assume you've aligned these sequences? I don't know GBlocks but I know trimal and will do the same thing (or do it better). Could you describe what you mean by 'removing gaps'? 1. Are you looking to trim the 5' and 3' ends? 2. Wanting to remove 'gap [indel] only columns'? 3. May be wanting to remove the any position which contains an indel, i.e. the entire site is removed? The last option is the best approach. $\endgroup$– M__ ♦Commented May 20, 2023 at 20:01
-
$\begingroup$ Did you look for the gblocks docs? You know it's >20 years old that tool! $\endgroup$– Chris_RandsCommented May 21, 2023 at 20:36
-
$\begingroup$ Hey @M__ ! I am still learning the procedures for doing phylogenetic analysis but I think I need to remove ambiguous sequences that hamper the accuracy of a phylogeny. You can check this paper which has cited Gblocks. [link] (doi.org/10.1016/j.syapm.2020.126166) If you feel trimal is better could you please explain how to use it? Thanks! $\endgroup$– K_081Commented May 22, 2023 at 5:44
-
$\begingroup$ @Chris_Rands, I did come across the Gblocks documentation but it did not work for me. I think it has not been updated in a while. Since my supervisor shared the paper which used Gblocks, I wanted to learn it. I agree that the latest tool would be better. $\endgroup$– K_081Commented May 22, 2023 at 5:50
Add a comment
|
1 Answer
$\begingroup$
$\endgroup$
2
Basically this is 16S data. If its proteins there's a different strategy - even if the data set needed is nucleotides (separate question). I would simply use trimal which is available on conda
trimal -in <inputfile> -out <outputfile> -nogaps
The input file is aligned fasta here and will output as aligned fasta.
Any site (alignment position) with an indel - will be removed. This is definitely preferred if performing an iqtree analysis, because of the way it handles indels. This will also clean up poorly aligned, alignments.
e.g.
---ATG
---A-G
--TATG
This will result in...
AG
AG
AG
Ultimately for alignments - there's no magic bullet - and manual editing and inspection is expected e.g. using seqotron. There are very good reasons for this concerning circular mathematics. However, modern aligners such as muscle5 are extremely good.
Trimal is available via conda
conda install -c bioconda trimal
In terms of using GBlock I'd need to know whether it was simply deleting indel sites, deleting indel sites at >50% frequency or whether it was performing a genetic diversity calculation per site. I don't think it is doing the latter. It is better to perform diversity saturation calculation explicitly rather than using an unknown automate because thats kinda important. If its deleting sites with excessive >50% indels - honestly its better to do a blockwise deletion rather than do something half-half. Within SSU there is a reason might be being done because there can be excessive indels at deep taxonomic levels with this data type. Again using iqtree just do it 'blockwise' which is the trimal output.
Note I've not read the original GBlock paper, but if you summarise its application I could consider it against trimal. However, I definitely prefer trimal considering the downstream application.
-
1$\begingroup$ Thank you again, @M__ . Will definitely go ahead with what you have suggested. Although, I must admit, being an absolute amateur in bioinformatics I don't follow some of the things you have so sincerely explained. Forgive me if I am being imprudent, do you have any way that I could learn all this under your guidance? I wish to be thorough.. thanks again! $\endgroup$– K_081Commented May 23, 2023 at 7:20
-
1$\begingroup$ If you need the theory you can always ask a separate question. The theory is very detailed, why certain things and done and otherwise not. $\endgroup$– M__ ♦Commented May 23, 2023 at 12:40