Feature: Check for hairpin artefact in linear bacterial chromosome ends

[Adamtaranto]

Some bacterial linear chromosomes terminate in a hairpin structure, whereby the top strand loops around to become its own complementary strand. This can result in reads that end in a spurious tandem inverted repeat (is actually just the initial strand being read backwards via the complementary strand).

See this figure from Autocycler: https://github.com/rrwick/Autocycler/wiki/Autocycler-dotplot#example-2-linear-plasmid

We could include a check and warning for this.

Alternatively, this could just be a recommended manual check.

Questions for @rrwick @dalofa

• What is the expected size range of these inverted repeats?
• Could they be the full chromosome length?
• Do we see hairpin reads in PacBio libs? I imagine the SMRT bell would only ligate to one end and this might confuse the CCS caller.
• Does this occur in all ONT library prep methods? I'd guess not in tagmentation libraries?

[dalofa]

I think checking whether the read (or must of the reads) chosen for extension contains a perfect hairpin onto itself (no spacer in between) would be a good way to check for it. I have attached an example from a hairpin end, where you can see that the reads start to map onto themselves from the soft-clipped part.

• So far I have observed from 8-9000 bases inverted repeat just from the reads directly. I think it is highly unlikely to observe them chromosome length, but definitely plasmid length. It seem to mostly relate to how far the nanopore continues to read. Also, there does not seem to be a drop in quality as the
• I have not yet seen any PacBio data from any of these strains, but will comment when/if I run into some.
• Regarding library prep: My expectations is that it can happen in both cases, but for now I know that it occurs in the ligation kit. I will comment when I have seen some of the rapid barcoding data.

An example of a hairpin end ONT ligation kit mapped to it:

[Adamtaranto]

Well that's interesting - if linear molecules always end in a hairpin and these always let ONT read back around onto the complementary strand, then the length distribution of the mirrored section should be the same as the distribution of read lengths overall. (or maybe with some depletion factor closer to the terminus due to fragment size selection like we were discussing).

In either case the reads have no reason to end at the hairpin point.

Which raises the question, why would the assembly contigs ever break exactly at the hairpin boundary? (as in your example)

Is that just the overhanging reads in the image? do you have the same dataset will all aligned reads visible at the contig end?

[dalofa]

This is from a manually curated assembly @rrwick made, so this is why the breakpoint is there. Typically, the issue will be in the initial assembly, where a "false" inverted repeat is created. So, I guess it makes more sense to look for the problem in the assembly rather than the reads. Warning-wise, I think it makes sense to implement, if you think it will be an issue that a user will take a curated assembly and then run teloclip on it.

I think you are right about the read distribution.

On the picture: There are more reads than what is easily shown, see below:

[Adamtaranto]

Ah that makes sense! So the IR is in the initial assembly and has been manually trimmed out.

Ok I will add a "check and warn" task to the todo list.

I'm undecided if this should happen in the extract module, or should be part of a suite of tests that is run before/after teloclip.

Proposed workflow:

• Extract last 10 Kbp of L/R end of each contig. (or whole contig if < 10Kbp) [Note this will miss IR if hairpin repeat is > 10Kbp. consider using whole chromosome]
• Self align with nucmer
• Check for inverted repeat that includes terminal section of contig
  • External IR begins no more that x bases from outer end,
  • Internal IR section is within x bases of external IR section (i.e. basically tandem).

[dalofa]

Follow-up here:
I have implemented a similar check in finishline using self-alignment from minimap2, which seems to work fine (although I am not sure the exact trimming position is always reached).

In some case the assemblies need to be trimmed with more than 10 kb (I have seen up to 30 kb). Without having though to much of it I think it might make sense to do num_of_bases_ext + 10-30 kb instead rather than just the terminal 10 kb.

[Adamtaranto]

Would the worst case be that the full linear chromosome (with hairpin closing off each end) is duplicated? I guess a single read could run around the full chromosome from a random breakpoint - behaves like a big ssDNA circle with complete tandem duplication.

 Linear chromosome with closed ends and a single-strand nick:

AAAAAAAAA/AAAAAAAA
[                  ]
 BBBBBBBBBBBBBBBBB

becomes:

AAAAAAAAA-BBBBBBBBBBBBBBBBB-AAAAAAAA

[dalofa]

In principle, yes.

I have never observed full chromosome duplication, but I have definitely seen plasmids that were full inverted repeats, which we suspected were assembly artifacts.