I think it will be necessary to filter out in some...
# phyloplace
a
I think it will be necessary to filter out in some way the output of hmmalign. Many of the results are mostly gaps and with the
--no-pre-mask
option it slows the overall performance. If I check the number of seqs in the output of hmmalign I can have something like this:
Copy code
file  format  type     num_seqs    sum_len  min_len  avg_len  max_len
-     FASTA   Protein   250,212  7,200,646        1     28.8      610
(take into account that the length is with gaps removed) If I filter out any sequence below 15 aa I get something like this:
Copy code
file  format  type     num_seqs    sum_len  min_len  avg_len  max_len
-     FASTA   Protein    66,256  6,846,372       15    103.3      610
It is quite a reduction. Could it be a step of the pipeline? Or else we should understand what is happening with hmmalign? Let me know what you think, thanks!
d
Filtering can certainly be done. I suppose one needs to do this after the profile alignment to the reference. Please make an issue and I'll look into this. Or if you want to get involved, make a PR! 🙂
I am, however, a little surprised that the sequences pass the HMM step. Are you using HMMs that are longer than the reference alignment?
a
I am using an HMM that is not directly equivalent to the alignment, perhaps this is the culprit, yes. But I am guessing that many users will be in this situation. The HMM I am using was used to find the sequences that afterwards are included in the alignment, but the alignment to do the tree was trimmed afterwards, and this may generate the difference? Given the normal user case, perhaps it would be necessary to include a comment on that. I will check it out and let you know!
d
If you by "trimmed" mean reduced to well-aligned positions I think your profiles are too wide since they appear to contain domains that are not used in the placement. I agree with you that this should be better documented in the pipeline. I also thought of adding the possibility to generate hmm files from the reference alignments which might be an easier way in some use cases. There is however a slight risk with that -- which might be relevant for your use case -- in that sometimes having only the target hmm in the file generates false positives of homologous proteins that are not covered by the tree. One sometimes needs baits in the form of hmms to attract away false positives. In your case, you seem to have the opposite problem: You're attracting too much. I'd recommend making profiles from the reference alignments (perhaps before you trimmed to well-aligned positions, but making sure they just contain the target domain) and use the current profiles as bait for false positives.