Build a core-genome phylogeny

Overview

This tutorial describes step by step the reconstruction of a core genome phylogeny for a given clade of bacteria.

Setup

1. clone eutils repository

   1 # download repo from metagenlab
   2 git clone https://github.com/metagenlab/eutils.git
   3 # move to eutils directory
   4 cd eutils
   5 # add path to pythonpath and path
   6 echo "export PATH=:$PWD"':$PATH' >> ~/.bashrc
   7 echo "export PYTHONPATH=:$PWD"':$PYTHONPATH' >> ~/.bashrc
   8 # reload bashrc
   9 source ~/.bashrc
   

2. install miniconda
3. create a dedicated environment and install orthofinder

   1# create a new environment
   2conda create -n orthofinder
   3source activate orthofinder
   4# install orthofinder
   5conda install orthofinder
   6# install biopython
   7conda install biopython
   

Get table of genome assemblies

The list of available genome assemblies can be downloaded with the script get_ncbi_genome_assemblies.py. Various options are available to filter the data in various ways:

• include only complete genomes
• include only representative or reference genomes (as classified on the NCBI)
• exclude incompletely classified genomes
• exclude genomes derived from metagenomes

If your taxon of interest includes too many genomes, it is also possible to get only one representative genome per species (or per genus/family/order). If you do so, one genome is selected among all assemblies available for each species. The genome is not selected randomly: priority is given to 1) complete genomes 2) reference genomes 3) representative genomes.

usage: get_ncbi_genome_assemblies.py [-h] [-t TAXON_ID] [-c] [-rf] [-rp] [-rs]
                                     [-gb] [-md] [-e] [-d] [-ra RANK]
                                     [-o OUTPUT]

optional arguments:
  -h, --help            show this help message and exit
  -t TAXON_ID, --taxon_id TAXON_ID
                        taxon_idd
  -c, --complete        Filter complete genomes (default: False)
  -rf, --reference      Filter reference genomes (default: False)
  -rp, --representative
                        Filter representative genomes (default: False)
  -rs, --refseq         Filter refseq assemblies (default: False)
  -gb, --genbank        Filter genbank assemblies (default: False)
  -md, --exclude_metagenome_derived
                        Exclude metagenome derived assemblies (default: False)
  -e, --exclude_unclassified
                        Exclude species labelled as "unclassified xxx" or as
                        "xxxx sp. yyy" (e.g. Blautia sp. TF11-31AT)
                        (default=False)
  -d, --complete_list   Report complete list of assemblies (no taxonomy
                        clustering) (default: False)
  -ra RANK, --rank RANK
                        Keep one prepresentative/clade of specified taxonomic
                        <rank>
  -o OUTPUT, --output OUTPUT
                        Output name (default: assemblies_<taxid>.tab )

Example: get one assembly per species for the Lachnospiraceae family

1. go on the NCBI taxonomy website. Search for Lachnospiraceae, then click to go to the main page of the family to find the unique taxonomic identifier of the family. It will look like this (the taxonomic ID is indicated with the red arrow):

1. execute the script using the family taxon ID. Since the Lachnospiraceae family encompass > 100 different species with multiple sequenced genome for each species, we will restrict the search to refseq genome, exclude metagenome derived assemblies, exclude incompletely classified genomes and only keep one genome per species:

1get_ncbi_genome_assemblies.py -t 186803 -rs -ra -e species -o lachnospiraceae_assemblies.tab

The search will take some time because the script perform a search for each species to identify the highest quality genome (complete genome whenever possible). In this particular case (in october 2018), this search reported a total of 706 genomes assemblies classified in 478 different taxon identifiers. After filtration, the number of genomes was reduced to 106 (keeping only 1 genome/species and excluding metagenome dervied assemblies and incompletely classified genomes).

The resulting table should look as follow:

rank_taxid | representative_species | taxid_representative_species | assembly_name | RefSeq | Genbank | contig_count | status | repres | length | clade | n_subtaxid_available | no_rank | superkingdom | phylum | class | order | family | genus |

-

28051 | Lachnospira multipara | 28051 | ASM42410v1 | GCA_000424105.1 | GCF_000424105.1 | 32 | Scaffold | representative genome | 2871860 | Lachnospira multipara | 8 | Terrabacteria group | Bacteria | Firmicutes | Clostridia | Clostridiales | Lachnospiraceae | Lachnospira |

-

1527 | Anaerocolumna aminovalerica | 1527 | IMG-taxon 2634166351 annotated assembly | GCA_900115365.1 | GCF_900115365.1 | 84 | Scaffold | representative genome | 4624179 | Anaerocolumna aminovalerica | 2 | Terrabacteria group | Bacteria | Firmicutes | Clostridia | Clostridiales | Lachnospiraceae | Anaerocolumna |