чoe*=:dZddlmZmZmZddlZddlmZddlZddl m Z ddl m Z ddlmZddlZddlmZddlmZddlZ ddlZd Zn #e$rd ZYnwxYwd Zd Zgd ZdZdZ d Z!dZ"dZ#gdZ$gdZ%dZ&dZ'dZ(dZ)ddZ*dZ+dZ,dZ-e.dkr e-dSdS)zo Create BIOM-format tables (http://biom-format.org) from Kraken output (http://ccb.jhu.edu/software/kraken/). )absolute_importdivisionprint_functionN) OrderedDict)datetime)open)dedent)TableTFShareef M. Dabdoubz"Copyright 2016, Shareef M. Dabdoub)r zAkshay ParopkarizSukirth Ganesanz Purnima KumarMITz'http://github.com/smdabdoub/kraken-biomzdabdoub.2@osu.eduz1.0.1) pct_reads clade_reads taxon_readsrankncbi_taxsci_name)DPCOFGSSScLdvrNdvrJddr)dtddzdd<dvrNdvrJddr)dtddzdd<fdtd|dzD}|dt|dzdD|d d rd |d ddz|d <|S) al Create a string representation of a taxonomic hierarchy (QIIME for now). :type tax_lvl: dict :param tax_lvl: Keyed on the entries in ranks :type end: int :param end: The end rank index (0-based indexing). If end == 3 then the returned string will contain K, P, C, and O. >>> tax_fmt({"K": "Bacteria", "P": "Firmicutes", "C": "Negativicutes", ... "O": "Selenomonadales", "F": "Veillonellaceae", "G": "Veillonella", ... "S": "Veillonella parvula", "SS": "Veillonella parvula DSM 2008"}, 4) 'k__Bacteria; p__Firmicutes; c__Negativicutes; o__Selenomonadales' rrNrcxg|]6}d||vr|nd7S)z{}__{}formatlower.0rtax_lvls `/mounts/lovelace/software/anaconda3/envs/kraken-biom/lib/python3.11/site-packages/kraken_biom.py ztax_fmt..CsO % % % ??17799ALLgajjb I I % % %c\g|])}d|*S)z{}__r)r#r$s r&r'ztax_fmt..Gs,BBBQ aggii((BBBr(rdk) startswithlenranksextend)r%endtaxs` r&tax_fmtr3,sN g~~ '>>gcl55gclCC>"3<GCL(9(9!(;(<(<=GCL w '>>gdm66ws|DD>#DM#gcl*;*;A*=*>*>?GDM % % % %FSUFm % % %CJJBB%Ab/BBBCCC  1v SVABBZA Jr(c t}t}i d}t|}t|}|D]o}|d} | tvrSt| |kr5t| } fdtd|D | tvr7|d | <t| }| tvr|t|dcxkr|krnt |d} t |d} | dks| dkrI|d|kr=t |||d<|d|kr | ||d<d| ||d<q||fS) a2 Parse a single output file from the kraken-report tool. Return a list of counts at each of the acceptable taxonomic levels, and a list of NCBI IDs and a formatted string representing their taxonomic hierarchies. :type kdata: str :param kdata: Contents of the kraken report file. rrc.i|]}||vr|ndS)rr"s r& z'parse_kraken_report..os+PPPqW '!**"PPPr(Nrrrr)rr/indexstripintr3) kdatamax_rankmin_ranktaxacountsr$ max_rank_idx min_rank_idxentryerankrrr%s @r&parse_kraken_reportrDTs ==D ]]FG A;;x((L;;x((L<<f ##%% E>>ekk%00144 E""APPPPeBQBiPPPG E>>":.4466GEN E""A E>>lekk%-.H.HXXXXLXXXXXeM233KeM233KQ;??uV}7P7P*1'1*=*=U:&'=H,,0;F5,--0;F5,- 4<r(ct}t}|D])}tj|s"td|tjtj|dd}t|d5} tj |td}d|Ddd} n4#t$r'} td | d} ~ wwxYw dddn #1swxYwYt| || \} } | | | ||<+||fS) zs Parse all kraken-report data files into sample counts dict and store global taxon id -> taxonomy data zERROR: File '{}' not found.rrrt ) fieldnames delimitercg|]}|Sr6r6)r#rBs r&r'z#process_samples..s0005000r(Nz ERROR: {}r<r=)rospisfile RuntimeErrorr splitextsplitrcsv DictReader field_namesOSErrorrDupdate) kraken_reports_fpr<r=r> sample_countskrep_fp sample_idkfkdrr;oescountsstaxas r&process_samplesr_s ==DMMM$++z'"" N<CCGLLMM ML7!3!3A!677: '4  ;B ;nRK/355500C0004 ; ; ;";#5#5b#9#9::: ;  ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;-UX6>@@@ E#* i  $ s6D0C  D C>"C99C>>DD D c vfdD}tj|t}fdD}dtt }t |tt|dttj |dS)a Create a BIOM table from sample counts and taxonomy metadata. :type sample_counts: dict :param sample_counts: A dictionary of dictionaries with the first level keyed on sample ID, and the second level keyed on taxon ID with counts as values. :type taxa: dict :param taxa: A mapping between the taxon IDs from sample_counts to the full representation of the taxonomy string. The values in this dict will be used as metadata in the BIOM table. :rtype: biom.Table :return: A BIOM table containing the per-sample taxon counts and full taxonomy identifiers as metadata for each taxon. c0g|]fdDS)cDg|]}|vrdn |S)rr6)r#sidrWtaxids r&r'z0create_biom_table...sF ( ( (}S111QQ}S7I%7P ( ( (r(r6)r#rdrWs @r&r'z%create_biom_table..sP # # # ( ( ( ( (& ( ( ( # # #r()dtypec$g|] }d|i S)taxonomyr6)r#rdr>s r&r'z%create_biom_table..s"<<datatax_metagen_strs`` r&create_biom_tablerxs  # # # #! # # #D 8D $ $ $D<<<EEeLLMMM b$'4 499W%%&&&''''''''''''''''''s')BB!$B!cd}tjt|tj}|ddd|ddt dd d |d d t dd d |dddd|dd|ddgdd |ddd|dd d!tt"|d#d$dd%| S)&Na3 Create BIOM-format tables (http://biom-format.org) from Kraken output (http://ccb.jhu.edu/software/kraken/). The program takes as input, one or more files output from the kraken-report tool. Each file is parsed and the counts for each OTU (operational taxonomic unit) are recorded, along with database ID (e.g. NCBI), and lineage. The extracted data are then stored in a BIOM table where each count is linked to the Sample and OTU it belongs to. Sample IDs are extracted from the input filenames (everything up to the '.'). OTUs are defined by the --max and --min arguments. By default these are set to Order and Species respectively. This means that counts assigned directly to an Order, Family, or Genus are recorded under the associated OTU ID, and counts assigned at or below the Species level are assigned to the OTU ID for the species. Setting a minimum rank below Species is not yet available. The BIOM format currently has two major versions. Version 1.0 uses the JSON (JavaScript Object Notation) format as a base. Version 2.x uses the HDF5 (Hierarchical Data Format v5) as a base. The output format can be specified with the --fmt option. Note that a tab-separated (tsv) output format is also available. The resulting file will not contain most of the metadata, but can be opened by spreadsheet programs. Version 2 of the BIOM format is used by default for output, but requires the Python library 'h5py'. If the library is not installed, kraken-biom will automatically switch to using version 1.0. Note that the output can optionally be compressed with gzip (--gzip) for version 1.0 and TSV files. Version 2 files are automatically compressed. Usage examples -------------- 1. Basic usage with default parameters: $ kraken-biom.py S1.txt S2.txt This produces a compressed BIOM 2.1 file: table.biom 2. BIOM v1.0 output: $ kraken-biom.py S1.txt S2.txt --fmt json Produces a BIOM 1.0 file: table.biom 3. Compressed TSV output: $ kraken-biom.py S1.txt S2.txt --fmt tsv --gzip -o table.tsv Produces a TSV file: table.tsv.gz 4. Change the max and min OTU levels to Class and Genus: $ kraken-biom.py S1.txt S2.txt --max C --min G Program arguments -----------------) descriptionformatter_classkraken_reports+z*Results files from the kraken-report tool.)nargshelpz--maxrr*zpAssigned reads will be recorded only if they are at or below max rank. Default: O.)defaultchoicesrz--minrzReads assigned at and below min rank will be recorded as being assigned to the min rank level. Default: S.z-oz --output_fpz table.biomaIPath to the BIOM-format file. By default, the table will be in the HDF5 BIOM 2.x format. Users can output to a different format using the --fmt option. The output can also be gzipped using the --gzip option. Default path is: ./table.biom)rrz--otu_fpaUCreate a file containing just the (NCBI) OTU IDs for use with a service such as phyloT (http://phylot.biobyte.de/) to generate a phylogenetic tree for use in downstream analysis such as UniFrac, iTol (itol.embl.de), or PhyloToAST (phylotoast.org).)rz--fmtry)ryr}rzVSet the output format of the BIOM table. Default is HDF5.z--gzip store_truezCompress the output BIOM table with gzip. HDF5 BIOM (v2.x) files are internally compressed by default, so this option is not needed when specifying --fmt hdf5.)actionrz --versionversionzkraken-biom version {}, {})rrz-vz --verbosezNPrints status messages during program execution.) argparseArgumentParsertwddRawDescriptionHelpFormatter add_argumentr/r rnro parse_args)descrparsers r&handle_program_optionsrs8 Et $e+3+OQQQF (IKKK eCRCjJKKK eCRCj;<<< m\?@@@   NOOO   7 7 70111 IJJJ  I188gNNPPP k,*+++     r(ct}|jdkr,ts%d|_d}tt |t |jt |jkr4d}tj | |j|jt|j |j|j\}}t||}t||j|j|j}|jrd t't)||jn@#t*$r3}d}tj | |Yd}~nd}~wwxYw|jrtd |d ||jd |jd | }tt |dSdS) Nryr}zo Library 'h5py' not found, unable to write BIOM 2.x (HDF5) files. Defaulting to BIOM 1.0 (JSON).z2ERROR: Max and Min ranks are out of order: {} < {}rKzERROR creating OTU file: {}rz BIOM-format table written to: {out_fp} Table contains {rows} rows (OTUs) and {cols} columns (Samples) and is {density:.1%} dense.rr)out_fprowscolsdensity)rr HAVE_H5PYprintrr/r8maxminsysexitr r_rrxrrrotu_fprrprNverboseshapeget_table_density)argsmsgrWr>rrre table_strs r&mainrds ! # #D x6)* d3ii {{48 DH 5 555B DHdh//000*$*=378378===M4 mT 2 2E t~tx C CF {% % 4::t{ 3 3 3 3 % % %3C HSZZ^^ $ $ $ $ $ $ $ $ % |  bii   '(.vV38;q>38;q>6;6M6M6O6O(.(Q(Q  d9oo  s"D:: E7)E22E7__main__)ryF)/__doc__ __future__rrrr collectionsrrQrrrrrros.pathpathrLrtextwrapr r biom.tabler numpyrlrr ImportError __author__ __copyright__ __credits__ __license__ro__maintainer__ __email__rnrSr/r3rDr_rxrrrr__name__r6r(r&rsA@@@@@@@@@###### ######"""""" ######KKKIIIII" 4 333  3%   /// 111%%%P333l@<<<:!!!!H'''&cccL)))X zDFFFFFsA AA