The goal with this analysis is to see how a given protein sequence is conserved across different taxonomy.

First we need to run blastp for our protein. To do this, we need to make sure blast command line tools are installed as described elsewhere. Personally I prefer to use conda to do this https://anaconda.org/bioconda/blast.

Now that this is set up, we need to blast our protein query and save the results of that search to file. Here, I am using the online database so that I don’t have to download the entire blast database locally, but this should really only be used for small searches as to not use to much of the NCBIs servers. Otherwise, they will cut off/limit your IP address eventually.

We can do this in bash as follows.

# See this link for discussion, it was very confusing at first to figure out how to get the taxonomy returned from the blastp call. - https://www.biostars.org/p/76551/


# Needed to download taxadb (converts numbers to taxonomy) and set the path variable to let blast know where to find it. 

# Downloaded taxdb via ftp@ ftp://ftp.ncbi.nlm.nih.gov/blast/db/taxdb.tar.gz
# Then transferred contents to /Users/charliebayne/blastdb/

#Setting path variable
export BLASTDB=$BLASTDB:/Users/charliebayne/blastdb/

#Looping through files
for i in *.fasta
  do blastp -db nr -query $i -remote -out $i.csv -outfmt "10 qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore staxids ssciname scomname sskingdom"
  done

#Combining all files
cat *.csv > blasp_hits_temp.csv
#Adding column names
echo -e "qseqid,sseqid,pident,length,mismatch,gapopen,qstart,qend,sstart,send,evalue,bitscore,staxids,ssciname,scomname,sskingdom" | cat - blasp_hits_temp.csv > blasp_hits.csv

#deleting temp
rm blasp_hits_temp

Note that this was somewhat confusing, as I needed to download the taxadb that allows for the conversions of NCBI taxonomy numbers to useful names and put it in the search path. Once I got that sorted outm it worked beautifully though.

Now we have a csv file containing all of the blastp hits. We can use that to plot the final taxonomic tree using the metacoder and taxize packages.

Now we can use a few R packages as described below to plot the final taxonomic tree.

library(metacoder)
library(taxize)
# Reading blastp data into R
trees_df = read_csv("blasp_hits.csv",col_names = TRUE)

# Getting higher classifications took code from https://www.r-bloggers.com/2018/11/using-taxa-and-metacoder-to-explore-taxonomy-of-vancouvers-trees/

species_df = as.data.frame(table(trees_df$ssciname))

species_classifications = classification(as.character(species_df$Var1), db = 'ncbi',rows = 1)
## ══  273 queries  ═════════════
## ✔  Found:  Alistipes+sp.
## ✔  Found:  Alistipes+sp.+Z76
## ✔  Found:  Alloprevotella+sp.
## ✔  Found:  Bacteroidaceae
## ✔  Found:  Bacteroidaceae+bacterium
## ✔  Found:  Bacteroidaceae+bacterium+14-104
## ✔  Found:  Bacteroidaceae+bacterium+HV4-6-C5C
## ✔  Found:  Bacteroidales+bacterium
## ✔  Found:  Bacteroidales+bacterium+Barb4
## ✔  Found:  Bacteroidales+bacterium+Barb6
## ✔  Found:  Bacteroidales+bacterium+Barb6XT
## ✔  Found:  Bacteroidales+bacterium+SW299
## ✔  Found:  Bacteroides
## ✔  Found:  Bacteroides+acidifaciens
## ✔  Found:  Bacteroides+bouchesdurhonensis
## ✔  Found:  Bacteroides+caccae
## ✔  Found:  Bacteroides+caecicola
## ✔  Found:  Bacteroides+caecigallinarum
## ✔  Found:  Bacteroides+caecimuris
## ✔  Found:  Bacteroides+cellulolyticus
## ✔  Found:  Bacteroides+cellulosilyticus
## ✔  Found:  Bacteroides+cellulosilyticus+DSM+14838
## ✔  Found:  Bacteroides+clarus
## ✔  Found:  Bacteroides+congonensis
## ✔  Found:  Bacteroides+cutis
## ✔  Found:  Bacteroides+eggerthii
## ✔  Found:  Bacteroides+faecalis
## ✔  Found:  Bacteroides+faecichinchillae
## ✔  Found:  Bacteroides+faecis
## ✔  Found:  Bacteroides+finegoldii
## ✔  Found:  Bacteroides+finegoldii+DSM+17565
## ✔  Found:  Bacteroides+fluxus
## ✔  Found:  Bacteroides+fragilis
## ✔  Found:  Bacteroides+fragilis+3_1_12
## ✔  Found:  Bacteroides+fragilis+638R
## ✔  Found:  Bacteroides+fragilis+CL03T00C08
## ✔  Found:  Bacteroides+fragilis+CL05T00C42
## ✔  Found:  Bacteroides+fragilis+CL07T00C01
## ✔  Found:  Bacteroides+fragilis+HMW+610
## ✔  Found:  Bacteroides+fragilis+HMW+615
## ✔  Found:  Bacteroides+fragilis+HMW+616
## ✔  Found:  Bacteroides+gallinaceum
## ✔  Found:  Bacteroides+gallinarum
## ✔  Found:  Bacteroides+graminisolvens
## ✔  Found:  Bacteroides+helcogenes
## ✔  Found:  Bacteroides+heparinolyticus
## ✔  Found:  Bacteroides+humanifaecis
## ✔  Found:  Bacteroides+ihuae
## ✔  Found:  Bacteroides+ilei
## ✔  Found:  Bacteroides+intestinalis
## ✔  Found:  Bacteroides+intestinalis+DSM+17393
## ✔  Found:  Bacteroides+luti
## ✔  Found:  Bacteroides+mediterraneensis
## ✔  Found:  Bacteroides+muris
## ✔  Found:  Bacteroides+ndongoniae
## ✔  Found:  Bacteroides+neonati
## ✔  Found:  Bacteroides+nordii
## ✔  Found:  Bacteroides+oleiciplenus
## ✔  Found:  Bacteroides+ovatus
## ✔  Found:  Bacteroides+ovatus+ATCC+8483
## ✔  Found:  Bacteroides+propionicigenes
## ✔  Found:  Bacteroides+pyogenes
## ✔  Found:  Bacteroides+salyersiae
## ✔  Found:  Bacteroides+salyersiae+CL02T12C01
## ✔  Found:  Bacteroides+sp
## ✔  Found:  Bacteroides+sp.
## ✔  Found:  Bacteroides+sp.+1_1_30
## ✔  Found:  Bacteroides+sp.+1001136B_160425_E2
## ✔  Found:  Bacteroides+sp.+1001302B_160321_D4
## ✔  Found:  Bacteroides+sp.+2_2_4
## ✔  Found:  Bacteroides+sp.+214
## ✔  Found:  Bacteroides+sp.+224
## ✔  Found:  Bacteroides+sp.+3_1_23
## ✔  Found:  Bacteroides+sp.+3_1_33FAA
## ✔  Found:  Bacteroides+sp.+3_1_40A
## ✔  Found:  Bacteroides+sp.+3_2_5
## ✔  Found:  Bacteroides+sp.+41_26
## ✔  Found:  Bacteroides+sp.+43_46
## ✔  Found:  Bacteroides+sp.+44_46
## ✔  Found:  Bacteroides+sp.+51
## ✔  Found:  Bacteroides+sp.+519
## ✔  Found:  Bacteroides+sp.+AF17-1
## ✔  Found:  Bacteroides+sp.+AF34-31BH
## ✔  Found:  Bacteroides+sp.+AM10-21B
## ✔  Found:  Bacteroides+sp.+AM16-13
## ✔  Found:  Bacteroides+sp.+AM16-15
## ✔  Found:  Bacteroides+sp.+AM16-24
## ✔  Found:  Bacteroides+sp.+AM23-18
## ✔  Found:  Bacteroides+sp.+AM37-9
## ✔  Found:  Bacteroides+sp.+An19
## ✔  Found:  Bacteroides+sp.+An269
## ✔  Found:  Bacteroides+sp.+An322
## ✔  Found:  Bacteroides+sp.+An51A
## ✔  Found:  Bacteroides+sp.+AR29
## ✔  Found:  Bacteroides+sp.+BFG-257
## ✔  Found:  Bacteroides+sp.+BFG-606
## ✔  Found:  Bacteroides+sp.+CBA7301
## ✔  Found:  Bacteroides+sp.+D2
## ✔  Found:  Bacteroides+sp.+D22
## ✔  Found:  Bacteroides+sp.+D55t1_190419_B8
## ✔  Found:  Bacteroides+sp.+ET225
## ✔  Found:  Bacteroides+sp.+ET336
## ✔  Found:  Bacteroides+sp.+GM023
## ✔  Found:  Bacteroides+sp.+HF-4919
## ✔  Found:  Bacteroides+sp.+K03
## ✔  Found:  Bacteroides+sp.+L10-4
## ✔  Found:  Bacteroides+sp.+M27
## ✔  Found:  Bacteroides+sp.+Marseille-P3684
## ✔  Found:  Bacteroides+sp.+Marseille-P8574
## ✔  Found:  Bacteroides+sp.+NM69_E16B
## ✔  Found:  Bacteroides+sp.+NSJ-21
## ✔  Found:  Bacteroides+sp.+NSJ-48
## ✔  Found:  Bacteroides+sp.+OF03-11BH
## ✔  Found:  Bacteroides+sp.+OF04-15BH
## ✔  Found:  Bacteroides+stercorirosoris
## ✔  Found:  Bacteroides+stercoris
## ✔  Found:  Bacteroides+thetaiotaomicron
## ✔  Found:  Bacteroides+thetaiotaomicron+dnLKV9
## ✔  Found:  Bacteroides+timonensis
## ✔  Found:  Bacteroides+togonis
## ✔  Found:  Bacteroides+uniformis
## ✔  Found:  Bacteroides+uniformis+CL03T00C23
## ✔  Found:  Bacteroides+uniformis+str.+3978+T3+i
## ✔  Found:  Bacteroides+xylanisolvens
## ✔  Found:  Bacteroides+zoogleoformans
## ✔  Found:  Bacteroidetes+bacterium
## ✔  Found:  Bacteroidetes+bacterium+HGW-Bacteroidetes-10
## ✔  Found:  Bacteroidetes+bacterium+oral+taxon+272+str.+F0290
## ✔  Found:  Barnesiella+sp.
## ✔  Found:  Candidatus+Avibacteroides+faecavium
## ✔  Found:  Candidatus+Bacteroides+intestinigallinarum
## ✔  Found:  Candidatus+Parabacteroides+faecavium
## ✔  Found:  Candidatus+Parabacteroides+intestinavium
## ✔  Found:  Candidatus+Paraprevotella+stercorigallinarum
## ✔  Found:  Candidatus+Phocaeicola+caecigallinarum
## ✔  Found:  Candidatus+Phocaeicola+excrementigallinarum
## ✔  Found:  Candidatus+Phocaeicola+excrementipullorum
## ✔  Found:  Candidatus+Phocaeicola+faecigallinarum
## ✔  Found:  Candidatus+Phocaeicola+faecipullorum
## ✔  Found:  Candidatus+Phocaeicola+gallinarum
## ✔  Found:  Candidatus+Phocaeicola+gallistercoris
## ✔  Found:  Candidatus+Phocaeicola+merdavium
## ✔  Found:  Candidatus+Phocaeicola+merdigallinarum
## ✔  Found:  Catenibacterium+mitsuokai
## ✔  Found:  Chitinophagaceae+bacterium
## ✔  Found:  Chitinophagales+bacterium
## ✔  Found:  Escherichia+coli
## ✔  Found:  Hydrotalea+sp.
## ✔  Found:  Jilunia+laotingensis
## ✔  Found:  Kiritimatiellae+bacterium
## ✔  Found:  Lacibacter+cauensis
## ✔  Found:  Lacibacter+luteus
## ✔  Found:  Lacibacter+sp.+S13-6-6
## ✔  Found:  Lentisphaerae+bacterium+GWF2_52_8
## ✔  Found:  Marinifilaceae+bacterium
## ✔  Found:  Mediterranea+massiliensis
## ✔  Found:  Mediterranea+sp.
## ✔  Found:  Mediterranea+sp.+ET5
## ✔  Found:  Muribaculaceae+bacterium
## ✔  Found:  Nubsella+zeaxanthinifaciens
## ✔  Found:  Odoribacter+sp.
## ✔  Found:  Odoribacter+sp.+43_10
## ✔  Found:  Odoribacter+splanchnicus
## ✔  Found:  Parabacteroides
## ✔  Found:  Parabacteroides+acidifaciens
## ✔  Found:  Parabacteroides+bouchesdurhonensis
## ✔  Found:  Parabacteroides+chinchillae
## ✔  Found:  Parabacteroides+distasonis
## ✔  Found:  Parabacteroides+faecis
## ✔  Found:  Parabacteroides+goldsteinii
## ✔  Found:  Parabacteroides+gordonii
## ✔  Found:  Parabacteroides+johnsonii
## ✔  Found:  Parabacteroides+massiliensis
## ✔  Found:  Parabacteroides+merdae
## ✔  Found:  Parabacteroides+provencensis
## ✔  Found:  Parabacteroides+sp.
## ✔  Found:  Parabacteroides+sp.+20_3
## ✔  Found:  Parabacteroides+sp.+52
## ✔  Found:  Parabacteroides+sp.+AD58
## ✔  Found:  Parabacteroides+sp.+AF17-28
## ✔  Found:  Parabacteroides+sp.+AF18-52
## ✔  Found:  Parabacteroides+sp.+AF19-14
## ✔  Found:  Parabacteroides+sp.+AF27-14
## ✔  Found:  Parabacteroides+sp.+AF39-10AC
## ✔  Found:  Parabacteroides+sp.+AF48-14
## ✔  Found:  Parabacteroides+sp.+AGMB00274
## ✔  Found:  Parabacteroides+sp.+AM08-6
## ✔  Found:  Parabacteroides+sp.+AM44-16
## ✔  Found:  Parabacteroides+sp.+An277
## ✔  Found:  Parabacteroides+sp.+CH2-D42-20
## ✔  Found:  Parabacteroides+sp.+GYB001
## ✔  Found:  Parabacteroides+sp.+HGS0025
## ✔  Found:  Parabacteroides+sp.+merdae-related_45_40
## ✔  Found:  Parabacteroides+sp.+N37
## ✔  Found:  Parabacteroides+sp.+TM07-1AC
## ✔  Found:  Parabacteroides+sp.+W1-Q-101
## ✔  Found:  Parabacteroides+sp.+ZJ-118
## ✔  Found:  Parabacteroides+timonensis
## ✔  Found:  Paracnuella+aquatica
## ✔  Found:  Paramuribaculum+intestinale
## ✔  Found:  Paraprevotella+sp.
## ✔  Found:  Paraprevotella+xylaniphila
## ✔  Found:  Pedobacter+sp.+ASV17
## ✔  Found:  Pedobacter+ureilyticus
## ✔  Found:  Phocaeicola
## ✔  Found:  Phocaeicola+abscessus
## ✔  Found:  Phocaeicola+barnesiae
## ✔  Found:  Phocaeicola+coprocola
## ✔  Found:  Phocaeicola+coprocola+DSM+17136
## ✔  Found:  Phocaeicola+coprophilus
## ✔  Found:  Phocaeicola+coprophilus+DSM+18228+=+JCM+13818
## ✔  Found:  Phocaeicola+dorei
## ✔  Found:  Phocaeicola+dorei+DSM+17855
## ✔  Found:  Phocaeicola+faecalis
## ✔  Found:  Phocaeicola+faecicola
## ✔  Found:  Phocaeicola+faecium
## ✔  Found:  Phocaeicola+fibrisolvens
## ✔  Found:  Phocaeicola+intestinalis
## ✔  Found:  Phocaeicola+massiliensis
## ✔  Found:  Phocaeicola+massiliensis+B84634+=+Timone+84634+=+DSM+17679+=+JCM+13223
## ✔  Found:  Phocaeicola+paurosaccharolyticus
## ✔  Found:  Phocaeicola+plebeius
## ✔  Found:  Phocaeicola+plebeius+DSM+17135
## ✔  Found:  Phocaeicola+salanitronis
## ✔  Found:  Phocaeicola+sartorii
## ✔  Found:  Phocaeicola+sp.
## ✔  Found:  Phocaeicola+vulgatus
## ✔  Found:  Prevotellaceae+bacterium
## ✔  Found:  Pseudocnuella+soli
## ✔  Found:  Rikenellaceae+bacterium
## ✔  Found:  Rufibacter+glacialis
## ✔  Found:  Rufibacter+sediminis
## ✔  Found:  Rufibacter+sp.+SYSU+D00344
## ✔  Found:  Rufibacter+sp.+SYSU+D00433
## ✔  Found:  Rufibacter+sp.+XAAS-G3-1
## ✔  Found:  Rufibacter+tibetensis
## ✔  Found:  Segetibacter+sp.+SYSU+D00508
## ✔  Found:  Tidjanibacter+sp.
## ✔  Found:  unclassified+Bacteroides
## ✔  Found:  unclassified+Parabacteroides
## ✔  Found:  uncultured+bacterium
## ══  Results  ═════════════════
## 
## • Total: 273 
## • Found: 241 
## • Not Found: 0
#Defining function to extract taxonomy from list
class_to_df <- function(x, ranking){
  if(!is.na(x)){ 
  return(x$name[x$rank == ranking])
    }
}

#Extracting taxonomy from different levels of taxonomy. 
species_df$superkingdom <- as.character(sapply(species_classifications, class_to_df, ranking="superkingdom"))

species_df$phylum <- as.character(sapply(species_classifications, class_to_df, ranking="phylum"))
species_df$order <- as.character(sapply(species_classifications, class_to_df, ranking="order"))
species_df$family <- as.character(sapply(species_classifications, class_to_df, ranking="family"))
species_df$genus <- as.character(sapply(species_classifications, class_to_df, ranking="genus"))
species_df$species <- as.character(sapply(species_classifications, class_to_df, ranking="species"))
species_df <- species_df[species_df$superkingdom != "NULL", ]

trees_df_class <- species_df
kept_columns <- c('Var1', 'phylum', 'order', 'family', 'genus', 'species')
trees_df_class <- trees_df_class[, kept_columns]

ranks <- c('phylum', 'order', 'family', 'genus', 'species')
obj <- parse_tax_data(trees_df_class, class_cols = ranks, named_by_rank = T)

# Constructing plot
tax_graph <- obj %>% 
  filter_taxa(nchar(taxon_names)!=0) %>%
  filter_taxa(taxon_ranks == "genus", supertaxa = TRUE) %>%
  heat_tree(node_label = taxon_names, 
            node_size = n_obs, 
            node_label_size = n_obs, 
            node_color = n_obs, 
            node_size_range = c(0.01, .05),  
            node_color_range = RColorBrewer::brewer.pal(5, "BuGn"), 
            node_label_size_range = c(0.02, 0.05),
            node_color_axis_label = "Number of Trees")

tax_graph

Figure 1. Here we can see that this protein we were interested in is primarily found within Bacteroides.