The Daffodil DNA Project: Results

This is the process of converting the signal from the MinIon sequencer into the DNA bases this represents. The instrument produces thousands of reads of varying length, up to 40-50kb in length. Basecalling was carried out using the Oxford Nanopore Guppy software (version 6.1.2) using the appropriate high accuracy model for the flowcells and sequencing kits used (dna_r9.4.1_450bps_hac).

Every base of sequence generated is assigned a quality score, which is a measure of how likely the base is to have been called correctly. These range from 0-40, using a logarithmic scale, where a score of 10 represents a 1 in 10 chance that the base call is incorrect, a score of 20 relating to a 1 in 100 chance, and 30 being 1 in 1000. Sequences are filtered using NanoFilt to remove sequences shorter than 300 bp and also those with a minimum average quality score below 10.

A large proportion of these samples does not consist of daffodil chloroplastic DNA, so to make the assembly process easier, contaminants are removed from the sequence data. Kraken 2.1.2 is also used to identify what the origin of the seqeunces in each sample. The sequence reads are mapped to a known daffodil chloroplast sequence using winnowmap 2.03, and those which do not have similarity with this reference sequence are discarded.

The individual sequence reads represent a part of the chloroplast genome, so a process called assembly is used to essentially find where all the reads overlap with each other to allow the sequence of the entire chloroplast genome to be determined. This produces sequences referred to as 'contigs', which contain contiguous stretches of sequence. In the ideal world we would have a single contig representing the entire genome, but in reality we typically see more than one. The assembly was carried out using Flye 2.9. There is typically far more data present in these samples than required for assembly, so the longest reads are used, with enough reads to cover the genome 50 times. Once assembled, the remaining reads are then mapped against the assembly to try to identify and remove errors in a process called polishing. 3 iterations of polishing were carried out in this case.

To determine the correct orientation and order of the contigs, they are aligned against an existing daffodil chloroplast genome, which also allows missing regions of the assembly to be identified. This process produces scaffolds, which are combine the contig sequences in the correct order and orientation, and fills any gaps in the sequence with 'N' bases, which represents an unknown base. Scaffolding was carried out using RagTag 2.1.0, with the MH706763.1 database sequence as a reference.

Since the chloroplast genome is a circular molecule, the linear representation produced by assembly may start at any point in the genome. To try to resolve this the sequence is aligned against the known chloroplast genome and the location of the start position identified. The sequence is cut at this point and the removed part added back to the end of the sequence

The position of the genes on the chloroplast genome was determined by an approach finding similiarity to those present in the MH706763 sequence. An initial annotation is carried out using the Rapid Annotation Transfer Tool, then once all sequences are available a more comprehensive annotation is carried out using GeSeq

The individual assembled chloroplast genomes are combined together and a multiple sequence alignment carried out using the Mafft aligner, with it's l-INS-i iterative pairwise local alignment method. Some assemblies have produced only partial assemblies, which if included in the phylogenetic analysis would greatly reduce the available usable sequence. Alignments are therefore carried out using all the sequences, and a separate alignment based on just chloroplast sequences of which >120kb has been assembled.

Phylogenetic analysis is carried out using RaxML with the GTRCAT substitution model. The analysis is bootstrapped to provide an estimate of the confidence of the placement of each sequence in the tree. The number of bootstrap iteriations is determed automatically using the 'autoMRE' method.

A plot of the resulting tree is produced using the R ggtree package.