KrakenTools is a suite of scripts to be used alongside the Kraken, KrakenUniq, Kraken 2, or Bracken programs. These scripts are designed to help Kraken users with downstream analysis of Kraken results.
For news and updates, refer to the github page: https://github.com/jenniferlu717/KrakenTools/
KrakenTools has been published on September 28, 2022 as part of a protocol paper for using the Kraken software suite. Please cite the following when using any KrakenTools script:
[Lu J, Rincon N, Wood D E, Breitwieser F P, Pockrandt C, Langmead B, Salzberg S L, Steinegger M. Metagenome analysis using the Kraken software suite. Nature Protocols, doi: 10.1038/s41596-022-00738-y (2022)] (https://www.nature.com/articles/s41596-022-00738-y)
Please also cite the relevant paper for usage of KrakenTools with any of the listed programs.
For issues with any of the above programs, please open a github issue on their respective github pages. This github repository is dedicated to only the scripts provided here.
- extract_kraken_reads.py
- combine_kreports.py
- kreport2krona.py
- kreport2mpa.py
- combine_mpa.py
- filter_bracken_out.py
- fix_unmapped.py
- make_ktaxonomy.py
- make_kreport.py
- alpha_diversity.py (see Diversity/README.md)
- beta_diversity.py (see Diversity/README.md)
Some of the scripts in this package will require installation of Biopython: https://biopython.org/. Otherwise, scripts should work with installation of python.
No installation required. All scripts are run on the command line as described.
Users can make scripts executable by running
chmod +x myscript.py
./myscript.py -h
This program extract reads classified at any user-specified taxonomy IDs. User must specify the Kraken output file, the sequence file(s), and at least one taxonomy ID. Additional options are specified below. As of April 19, 2021, this script is compatible with KrakenUniq/Kraken2Uniq reports.
python extract_kraken_reads.py
-k, --kraken MYFILE.KRAKEN.............Kraken output file-s, -s1, -1, -U SEQUENCE.FILE..........FASTA/FASTQ sequence file (may be gzipped)-s2, -2 SEQUENCE2.FILE.................FASTA/FASTQ sequence file (for paired reads, may be gzipped)-o, --output2 OUTPUT.FASTA.............output FASTA/Q file with extracted seqs-t, --taxid TID TID2 etc...............list of taxonomy IDs to extract (separated by spaces)
Optional:
-o2, --output2 OUTPUT.FASTA.............second output FASTA/Q file with extracted seqs (for paired reads)--fastq-output..........................Instead of producing FASTA files, print FASTQ files (requires FASTQ input)--exclude...............................Instead of finding reads matching specified taxids, finds reads NOT matching specified taxids.-r, --report MYFILE.KREPORT.............Kraken report file (required if specifying --include-children or --include-parents)--include-children......................include reads classified at more specific levels than specified taxonomy ID levels.--include-parents.......................include reads classified at all taxonomy levels between root and the specified taxonomy ID levels.--max #.................................maximum number of reads to save.--append................................if output file exists, appends reads--noappend..............................[default] rewrites existing output file
Input sequence files must be either FASTQ or FASTA files. Input files can be gzipped or not. The program will automatically detect whether the file is gzipped and whether it is FASTQ or FASTA formatted based on the first character in the file (">" for FASTA, "@" for FASTQ)
Users that ran Kraken using paired reads should input both read files into extract_kraken_reads.py as follows:
extract_kraken_reads.py -k myfile.kraken -s1 read1.fq -s2 reads2.fq
Given paired reads, the script requires users to provide two output file names to contain extracted reads:
extract_kraken_reads.py -k myfile.kraken -s1 read1.fq -s2 reads2.fq -o extracted1.fq -o2 extracted2.fq
The delimiter (--delimiter or -d) option has been removed.
`extract_kraken_reads.py -k myfile.kraken ... -o reads_S1.fa -o2 reads_s2.fa
By default, reads classified at specified taxonomy IDs will be extracted (and any taxids selected using --include-parents/--include-children. However, specifying --exclude will cause the reads NOT classified at any specified taxonomy IDs.
For example:
extract_kraken_reads.py -k myfile.kraken ... --taxid 9606 --exclude==> extract all reads NOT classified as Human (taxid 9606).extract_kraken_reads.py -k myfile.kraken ... --taxid 2 --exclude --include-children==> extract all reads NOT classified as Bacteria (taxid 2) or any classification in the Bacteria subtree.extract_kraken_reads.py -k myfile.kraken ... --taxid 9606 --exclude --include-parents==> extract all reads NOT classified as Human or any classification in the direct ancestry of Human (e.g. will exclude reads classified at the Primate, Chordata, or Eukaryota levels).
By default, only reads classified exactly at the specified taxonomy IDs will be extracted. Options --include-children and --include parents can be used to extract reads classified within the same lineage as a specified taxonomy ID. For example, given a Kraken report containing the following:
[%] [reads] [lreads][lvl] [tid] [name]
100 1000 0 R 1 root
100 1000 0 R1 131567 cellular organisms
100 1000 50 D 2 Bacteria
0.95 950 0 P 1224 Proteobacteria
0.95 950 0 C 1236 Gammaproteobacteria
0.95 950 0 O 91347 Enterobacterales
0.95 950 0 F 543 Enterobacteriaceae
0.95 950 0 G 561 Escherichia
0.95 950 850 S 562 Escherichia coli
0.05 50 50 S1 498388 Escherichia coli C
0.05 50 50 S1 316401 Escherichia coli ETEC
extract_kraken_reads.py [options] -t 562==> 850 reads classified as E. coli will be extractedextract_kraken_reads.py [options] -t 562 --include-parents==> 900 reads classified as E. coli or Bacteria will be extractedextract_kraken_reads.py [options] -t 562 --include-children==> 950 reads classified as E. coli, E. coli C, or E. coli ETEC will be extractedextract_kraken_reads.py [options] -t 498388==> 50 reads classified as E. coli C will be extractedextract_kraken_reads.py [options] -t 498388 --include-parents==> 950 reads classified as E. coli C, E. coli, or Bacteria will be extractedextract_kraken_reads.py [options] -t 1 --include-children==> All classified reads will be extracted
This script combines multiple Kraken reports into a combined report file.
python complete_kreports.py
-r 1.KREPORT 2.KREPORT........................Kraken-style reports to combine-o COMBINED.KREPORT...........................Output file
Optional:
--display-headers..............................include headers describing the samples and columns [all headers start with #]--no-headers...................................do not include headers in output--sample-names.................................give abbreviated names for each sample [default: S1, S2, ... etc]--only-combined................................output uses exact same columns as a single Kraken-style report file. Only total numbers for read counts and percentages will be used. Reads from individual reports will not be included.
Percentage is only reported for the summed read counts, not for each individual sample.
The output file therefore contains the following tab-delimited columns:
perc............percentage of total reads rooted at this cladetot_all ........total reads rooted at this clade (including reads at more specific clades)tot_lvl.........total reads at this clade (not including reads at more specific clades)1_all...........reads from Sample 1 rooted at this clade1_lvl...........reads from Sample 1 at this clade2_all...........""2_lvl...........""- etc..
lvl_type........Clade level type (R, D, P, C, O, F, G, S....)taxid...........taxonomy ID of this cladename............name of this clade
This program takes a Kraken report file and prints out a krona-compatible TEXT file
python kreport2krona.py
-r/--report MYFILE.KREPORT........Kraken report file-o/--output MYFILE.KRONA..........Output Krona text file
Optional:
--no-intermediate-ranks...........[default]only output standard levels [D,P,C,O,F,G,S]--intermediate-ranks..............include non-standard levels
kraken2 --db KRAKEN2DB --threads THREADNUM --report MYSAMPLE.KREPORT \
--paired SAMPLE_1.FASTA SAMPLE_2.FASTA > MYSAMPLE.KRAKEN2
python kreport2krona.py -r MYSAMPLE.KREPORT -o MYSAMPLE.krona
ktImportText MYSAMPLE.krona -o MYSAMPLE.krona.html
Krona information: see https://github.com/marbl/Krona.
--no-intermediate-ranks
6298 Unclassified
8 k__Bacteria
4 k__Bacteria p_Proteobacteria
6 k__Bacteria p_Proteobacteria c__Gammaproteobacteria
...
--intermediate-ranks
6298 Unclassified
79 x__root
0 x__root x__cellular_organisms
8 x__root x__cellular organisms k__Bacteria
4 x__root x__cellular organisms k__Bacteria p__Proteobacteria
6 x__root x__cellular organisms k__Bacteria p__Proteobacteria c__Gammaproteobacteria
....
This program takes a Kraken report file and prints out a mpa (MetaPhlAn) -style TEXT file
python kreport2mpa.py
-r/--report MYFILE.KREPORT........Kraken report file-o/--output MYFILE.MPA.TXT........Output MPA-STYLE text file
Optional:
--display-header..................display header line (#Classification, MYFILE.KREPORT) [default: no header]--no-intermediate-ranks...........[default] only output standard levels [D,P,C,O,F,G,S]--intermediate-ranks..............include non-standard levels--read-count......................[default] use read count for output--percentages.....................use percentage of total reads for output--remove-spaces...................[default] replace spaces with underscores in each taxon--keep-spaces.....................keep spaces in each taxon
kraken2 --db KRAKEN2DB --threads THREADNUM --report MYSAMPLE.KREPORT \
--paired SAMPLE_1.FASTA SAMPLE_2.FASTA > MYSAMPLE.KRAKEN2
python kreport2mpa.py -r MYSAMPLE.KREPORT -o MYSAMPLE.MPA.TXT
The output will contain one tab character inbetween the classification and the read count.
--no-intermediate-ranks/--read-count
#Classification SAMPLE.KREPORT
k__Bacteria 36569
k__Bacteria|p__Proteobacteria 21001
k__Bacteria|p__Proteobacteria|c__Gammaproteobacteria 11648
...
--intermediate-ranks/--read-count
#Classification SAMPLE.KREPORT
x__cellular_organisms 38462
x__cellular_organisms|k__Bacteria 36569
x__cellular_organisms|k__Bacteria|p__Proteobacteria 21001
...
This program combines multiple outputs from kreport2mpa.py. Files to be combined must have been generated using the same kreport2mpa.py options.
Important:
- Input files to combine_mpa.py cannot be a mix of intermediate/no intermediate rank outputs.
- Input files should be generated using the same Kraken database.
- Input files cannot be a mix of read counts/percentage kreport2mpa.py outputs. combine_mpa.py will not test the input files prior to combining.
If no header is in a given sample file, the program will number the files "Sample #1", "Sample #2", etc.
python combine_mpa.py
-i/--input MYFILE1.MPA MYFILE2.MPA.......Multiple MPA-STYLE text files (separated by spaces)-o/--output MYFILE.COMBINED.MPA..........Output MPA-STYLE text file
#Classification Sample #1 Sample #2
k__Bacteria 36569 20034
k__Bacteria|p__Proteobacteria 21001 18023
k__Bacteria|p__Proteobacteria|c__Gammaproteobacteria 11648 15000
This program takes the output file of a Bracken report and filters the desired taxonomy IDs.
python filter_bracken_out.py
-i/--input MYFILE.BRACKEN..........Bracken output file-o/--output MYFILE.BRACKEN_NEW.....Bracken-style output file with filtered taxids--include TID TID2.................taxonomy IDs to include in output file [space-delimited]--exclude TID TID2.................taxonomy IDs to exclude in output file [space-delimited]
User should specify either taxonomy IDs with --include or --exclude. If
both are specified, taxonomy IDs should not be in both lists and only
taxonomies to include will be evaluated.
When specifying the --include flag, only lines for the included taxonomy IDs will be extracted to the filtered output file. The percentages in the filtered file will be re-calculated so the total percentage in the output file will sum to 100%.
When specifying the --exclude flag alone, all lines in the Bracken file will be preserved EXCEPT for the lines matching taxonomy IDs provided.
This program can be useful for isolating a subset of species to better understand the distribution of those particular species in the sample.
For example:
python filter_bracken_out.py [options] --include 1764 1769 1773 1781 39689will allow users to get the relative percentages of Mycobacterium avium, marinum, tuberculosis, leprae, and gallinarum in their samples.
In other cases, users may want to focus on the distribution of all species that are NOT the host species in a given sample. This program can then recalculate percentage distributions for species when excluding reads for the host.
For example, given this output:
name tax_id tax_lvl kraken.... added... new.... fraction...
Homo sapiens 9606 S ... .... 999000 0.999000
Streptococcus pyogenes 1314 S ... .... 10 0.000001
Streptococcus agalactiae 1311 S ... .... 5 0.000000
Streptococcus pneumoniae 1313 S ... .... 3 0.000000
Bordetella pertussis 520 S ... .... 20 0.000002
...
Users may not be interested in the 999,000 reads that are host DNA, but
would rather know the percentage of non-host reads for each of the non-host
species. Using python filter_bracken_out.py [options] --exclude 9606
allows better resolution of the non-host species, allowing each of the
fraction of reads to be recalculated out of 1,000 instead of 1,000,000
reads in the above example. The output would then be:
name tax_id tax_lvl kraken.... added... new.... fraction...
Streptococcus pyogenes 1314 S ... .... 10 0.01000
Streptococcus agalactiae 1311 S ... .... 5 0.05000
Streptococcus pneumoniae 1313 S ... .... 3 0.03000
Bordetella pertussis 520 S ... .... 200 0.20000
...
When building a Kraken database, an "unmapped.txt" file may be generated if a taxonomy for a given sequence is not found. This script can search through any accession2taxid files provided and the unmapped.txt file and generate a seqid2taxid.map file to be appended to the one already generated.
python fix_unmapped.py
-i/--input unmapped.txt...........Any file containing accession IDs to map--accession2taxid REF_FILES.......Any tab-delimited file with 4 columns, accessions = column 1, taxonomy IDs = column 3-o/--output OUT_FILE..............Output tab-delimited file with 2 columns: accessions and taxids
Optional:
-r/--remaining....................file containing any unmapped accession IDs after search [default:still_unmapped.txt]
rm *.k2d
mv seqid2taxid.map seqid2taxid_1.map
python fix_unmapped -i unmapped.txt --accession2taxid taxonomy/*accession2taxid -o seqid2taxid_temp.map
cat seqid2taxid_1.map seqid2taxid_temp.map
kraken2-build --build --db . --threads 4
For future KrakenTools scripts, this program generates a single text file that contains all of the taxonomy information required. This program is intended to generate a single text taxonomy file for any Kraken 1, Kraken 2, or KrakenUniq database.
Important: The output of this program does not replace any Kraken database file (do not replace your taxo.k2d or .db files).
python make_ktaxonomy.py
--nodes taxonomy/nodes.dmp...........nodes.dmp file in Kraken DB taxonomy/ folder--names taxonomy/names.dmp...........names.dmp file in Kraken DB taxonomy/ folder--seqid2taxid seqid2taxid.map........seqid2taxid.map file generated by kraken-build/kraken2-build/krakenuniq-build when building the database. This is a 2-column tab-delimited file containing sequence IDs and taxonomy IDs.-o/--output OUT_FILE.................Output text file. More details below
The program will inform users if a taxonomy ID is listed in the seqid2taxid.map
file but not in either the nodes.dmp or the names.dmp files.
The output file is similar to the nodes.dmp/names.dmp file format, but not identical.
Each of the following columns is separated by a tab-vertical line-tab (e.g. \t|\t).
- taxonomy ID
- parent taxonomy ID
- rank type (R = root, D = domain/superkingdom, P = phylum, etc.)
- level number (distance from root)
- name
For ranks outside of the traditional taxonomy ranks (R, D, P, C, O, F, G, S),
the rank type will be assigned based on the closest parent, with a number to specify
distance from that parent. For example, the strains will be labeled with S1 while
ranks inbetween Genus and Species will be labeled with G1, G2, etc.
Currently, names for each node are selected based on the first name listed in
the names.dmp file or the name designated as scientific name.
scientific names will be preferred over all others.
- taxonomy/nodes.dmp
- taxonomy/names.dmp
- seqid2taxid.map
This program will generate a kraken-style report file from the kraken output file. Currently, this only generates reports for Kraken 1 or Kraken 2. This program does not currently work for KrakenUniq output files (to be completed in a future project).
This program requires that users first generate the taxonomy file created by make_ktaxonomy.py.
python make_kreport.py
-i/-k/--input KRAKEN_FILE........default Kraken output file (5 tab-delimited columns, taxid in third column)-t/--taxonomy TAXONOMY_FILE......output from make_ktaxonomy.py-o/--output REPORT_FILE..........output Kraken report file (6 tab-delimited columns)
Optional
--use-read-len...................make report using summed read lengths instead of read counts
Given a Kraken 2 database KRAKENDB/ and sample file EXAMPLE_READS.fq,
the following commands can be used to generate a Kraken report file
with this script.
python make_ktaxonomy.py --nodes KRAKENDB/taxonomy/nodes.dmp --names KRAKENDB/taxonomy/names.dmp --seqid2taxid KRAKENDB/seqid2taxid.map -o KRAKENDB/mydb_taxonomy.txt
kraken2 --db KRAKENDB --threads 4 EXAMPLE_READS.fq > EXAMPLE.kraken2
python make_kreport.py -i EXAMPLE.kraken2 -t KRAKENDB/mydb_taxonomy.txt -o EXAMPLE.kreport2
By default, the output Kraken report will list read counts for each taxonomy ID. However,
if all read lengths are not the same, users can add the --use-read-len option, which will
result in reporting summed read lengths for each taxon.
The output format for kreport.py is identical to the format generated by
kraken-report or the --report switch with kraken2. The output
file contains 6 tab-delimited columns as follows:
- Percentage of total reads
- Reads classified within sub-tree
- Reads classified at this specific node (reads cannot be more specifically classified)
- Level type (R = root, K = kingdom, P = phylum, etc)
- Taxonomy ID
- Name (preceeded by spaces to indicate distance from root)
Jennifer Lu [email protected] [email protected] Page Updated: 2020/05/10
