Answers.Body: <p>I need to do a more thorough testing to determine just how effective this answer may be, but I think I've put together a workflow that is close enough to an appropriate answer that I'd post it here.</p>
  <p>In my research I came across <a href="https://almob.biomedcentral.com/articles/10.1186/s13015-018-0130-7" rel="nofollow noreferrer">gappy</a>, (GitHub <a href="https://github.com/alexdonath/gappy2" rel="nofollow noreferrer">here</a>), a piece of software that identifies potentially relevant gaps in a phylogenetic analysis and returns a PHYLIP-formatted file that can be run through most tree-building programs. Gappy accounts for indels, but not substitutions, so the output would need to be combined with another analysis to generate the tree I want.</p>
  <p>Fortunately, IQTree allows for a partitioned analysis of mixed data (see <a href="http://www.iqtree.org/doc/Advanced-Tutorial" rel="nofollow noreferrer">here</a>), so by feeding the original alignment AND the gappy output through IQTree you should be able to generate a tree that accounts for both substitutions and indels.</p>
  <p>Again, more thorough testing is required, but I tested this possibility using an MSA example taken from <a href="https://www.ebi.ac.uk/seqdb/confluence/display/JDSAT/Multiple+Sequence+Alignment+Tool+Output+Examples" rel="nofollow noreferrer">this website</a>. Saving the FASTA-formatted example as <code>alignment.fasta</code>, and running that file through gappy2, you get an output file called <code>GAPPY2_alignment.fasta_2-inf.phy</code>. Thus if you create a NEXUS-formatted file like the following:</p>
  <pre><code>#nexus
  begin sets;
      charset part1 = alignment.fasta: *;
      charset part2 = GAPPY2_alignment.fasta_2-inf.phy: *;
  end;
  </code></pre>
  <p>and run it through IQTree using the command <code>iqtree -p example.nex</code>, you get the following tree:</p>
  <pre><code>
  +---ENA|CAA23748|CAA23748.1
  |
  +-------------ENA|CAA24095|CAA24095.1
  |
  |              +-------------------------------------------ENA|BAA20512|BAA20512.1
  +--------------|
                 +-----------------------ENA|CAA28435|CAA28435.1
  </code></pre>
  <p>Now, disclaimer: this is the exact same tree you get when you run <code>alignment.fasta</code> through IQTree, hence why more thorough testing is required. The sample size for this example is quite small, and it's therefore not surprising that accounting for indels does not change the outcome in this scenario. However, I think this serves as a solid baseline workflow, where future updates will be more tweaks than complete overhauls, if that makes sense.</p>
  <p>Let me know if you see any glaring flaws in my logic, but I think I'm going to close this question. I may provide edits to this answer in the future if significant tweaks/overhauls wind up being done.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-05-11T14:19:58.827Z
Answers.LastActivityDate: 2022-05-11T14:19:58.827Z
Answers.UserDisplayName: CCranney
Body: <p>I'm part of a nanopore sequencing experiment that will sequence several generations of viruses. The intent is to perform directed evolution by putting selective pressure on these viruses and tracking the various mutations that occur. See <a href="https://www.cell.com/cell/fulltext/S0092-8674(19)30622-1?_returnURL=https%3A%2F%2Flinkinghub.elsevier.com%2Fretrieve%2Fpii%2FS0092867419306221%3Fshowall%3Dtrue" rel="nofollow noreferrer">here</a> for the paper this methodology is based on.</p>
  <p>I've been working on a software package that can compile consensus sequences from longread nanopore sequencing (see <a href="https://github.com/JGEnglishLab/longread_umi_python" rel="nofollow noreferrer">here</a>). However, there's a wrinkle with generating phylogenetic trees from this type of data. It appears that most tree building software ignores indels. This makes sense, as gaps could easily indicate a lack of coverage from illumina sequencing methods or similar issues rather than betray the existence of an insertion/deletion. However, the nature of our sequencing method is such that a gaps really should be evaluated as mutational changes that we want to keep track of in a phylogenetic tree.</p>
  <p>We've explored using NextStrain to generate our trees and visuals, but the tree building algorithms they employ too suffer from this gap assumption (see the issue I raised <a href="https://github.com/nextstrain/auspice/issues/1500" rel="nofollow noreferrer">here</a>). These include IQTree, RAxML and PhyML. I have not had much luck finding an algorithm or software that treats gaps as indels when building a phylogenetic tree. Does anyone know of such an algorithm or package I could employ? I figure I could take a stab at writing my own algorithm, but consider that a Plan B resort. No sense re-inventing the wheel and all that.</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2022-04-18T03:04:21.460Z
Comments.Text: Speaking personally it will take me a good half day to understand the VEGAS system (Cell paper), which - whilst I get the gist: using viruses as replication vectors to accelerate adaptation, in particularly I think it leverages copy number thresholds of capsid proteins; I don't understand the nuance of the wet-lab system.
Comments.UserDisplayName: M__
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2022-04-16T00:41:36.237Z
Comments.Text: Various distance measures can do this given an MSA. See for example https://academic.oup.com/bioinformatics/article/28/4/495/212883. As long as you have a transition matrix that deals with gaps it should work fine. You can build trees with e.g. neighbor joining from distance matrices. References on ML methods for dealing with indels can be found in [this paper](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3299439/).
Comments.UserDisplayName: Maximilian Press
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2022-04-16T03:26:07.570Z
Comments.Text: What the OP is saying is correct wft mainstream phylogenetics: the reason is they model point mutations - not alignment reliability. There are numerous methods (> 2 contrasting methods) to do what the OP seeks, its quite easy. There's no need to write their own algorithm and would be discourage (account for lots of parameters). Generic NJ is okish (its pairwise), its not quite what the OP wants. The virus cited in the Cell paper is Sindbis virus (alphavirus) ... very interesting. I think I'm probably going to sit out on this question, but I definitely look forward to seeing the paper. Upvoted.
Comments.UserDisplayName: M__
Tags: <phylogenetics><nanopore><long-reads><indel>
Title: Tree Building Algorithm that treats gaps as deletions
ThreadText: # Tree Building Algorithm that treats gaps as deletions
  
  (Asked by: username_0 on 2022-04-15)
  
  I'm part of a nanopore sequencing experiment that will sequence several generations of viruses. The intent is to perform directed evolution by putting selective pressure on these viruses and tracking the various mutations that occur. See [here](https://www.cell.com/cell/fulltext/S0092-8674(19)30622-1?_returnURL=https%3A%2F%2Flinkinghub.elsevier.com%2Fretrieve%2Fpii%2FS0092867419306221%3Fshowall%3Dtrue) for the paper this methodology is based on.
  
  I've been working on a software package that can compile consensus sequences from longread nanopore sequencing (see [here](https://github.com/JGEnglishLab/longread_umi_python)). However, there's a wrinkle with generating phylogenetic trees from this type of data. It appears that most tree building software ignores indels. This makes sense, as gaps could easily indicate a lack of coverage from illumina sequencing methods or similar issues rather than betray the existence of an insertion/deletion. However, the nature of our sequencing method is such that a gaps really should be evaluated as mutational changes that we want to keep track of in a phylogenetic tree.
  
  We've explored using NextStrain to generate our trees and visuals, but the tree building algorithms they employ too suffer from this gap assumption (see the issue I raised [here](https://github.com/nextstrain/auspice/issues/1500)). These include IQTree, RAxML and PhyML. I have not had much luck finding an algorithm or software that treats gaps as indels when building a phylogenetic tree. Does anyone know of such an algorithm or package I could employ? I figure I could take a stab at writing my own algorithm, but consider that a Plan B resort. No sense re-inventing the wheel and all that.
  
  **Tags:** phylogenetics, nanopore, long-reads, indel
  
  > **Comments:**
  > **username_2**: Various distance measures can do this given an MSA. See for example https://academic.oup.com/bioinformatics/article/28/4/495/212883. As long as you have a transition matrix that deals with gaps it should work fine. You can build trees with e.g. neighbor joining from distance matrices. References on ML methods for dealing with indels can be found in [this paper](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3299439/).  
  > **username_1**: What the OP is saying is correct wft mainstream phylogenetics: the reason is they model point mutations - not alignment reliability. There are numerous methods (> 2 contrasting methods) to do what the OP seeks, its quite easy. There's no need to write their own algorithm and would be discourage (account for lots of parameters). Generic NJ is okish (its pairwise), its not quite what the OP wants. The virus cited in the Cell paper is Sindbis virus (alphavirus) ... very interesting. I think I'm probably going to sit out on this question, but I definitely look forward to seeing the paper. Upvoted.  
  > **username_1**: Speaking personally it will take me a good half day to understand the VEGAS system (Cell paper), which - whilst I get the gist: using viruses as replication vectors to accelerate adaptation, in particularly I think it leverages copy number thresholds of capsid proteins; I don't understand the nuance of the wet-lab system.  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_0 on 2022-05-11. Score: 1 | Accepted answer
  
  I need to do a more thorough testing to determine just how effective this answer may be, but I think I've put together a workflow that is close enough to an appropriate answer that I'd post it here.
  
  In my research I came across [gappy](https://almob.biomedcentral.com/articles/10.1186/s13015-018-0130-7), (GitHub [here](https://github.com/alexdonath/gappy2)), a piece of software that identifies potentially relevant gaps in a phylogenetic analysis and returns a PHYLIP-formatted file that can be run through most tree-building programs. Gappy accounts for indels, but not substitutions, so the output would need to be combined with another analysis to generate the tree I want.
  
  Fortunately, IQTree allows for a partitioned analysis of mixed data (see [here](http://www.iqtree.org/doc/Advanced-Tutorial)), so by feeding the original alignment AND the gappy output through IQTree you should be able to generate a tree that accounts for both substitutions and indels.
  
  Again, more thorough testing is required, but I tested this possibility using an MSA example taken from [this website](https://www.ebi.ac.uk/seqdb/confluence/display/JDSAT/Multiple+Sequence+Alignment+Tool+Output+Examples). Saving the FASTA-formatted example as `alignment.fasta`, and running that file through gappy2, you get an output file called `GAPPY2_alignment.fasta_2-inf.phy`. Thus if you create a NEXUS-formatted file like the following:
  
  ```
  #nexus
  begin sets;
      charset part1 = alignment.fasta: *;
      charset part2 = GAPPY2_alignment.fasta_2-inf.phy: *;
  end;
  
  ```
  
  and run it through IQTree using the command `iqtree -p example.nex`, you get the following tree:
  
  ```
  
  +---ENA|CAA23748|CAA23748.1
  |
  +-------------ENA|CAA24095|CAA24095.1
  |
  |              +-------------------------------------------ENA|BAA20512|BAA20512.1
  +--------------|
                 +-----------------------ENA|CAA28435|CAA28435.1
  
  ```
  
  Now, disclaimer: this is the exact same tree you get when you run `alignment.fasta` through IQTree, hence why more thorough testing is required. The sample size for this example is quite small, and it's therefore not surprising that accounting for indels does not change the outcome in this scenario. However, I think this serves as a solid baseline workflow, where future updates will be more tweaks than complete overhauls, if that makes sense.
  
  Let me know if you see any glaring flaws in my logic, but I think I'm going to close this question. I may provide edits to this answer in the future if significant tweaks/overhauls wind up being done.
  
  
  ---
  
  