#!/usr/bin/env python3

# supress python warnings
import warnings
#from flask.exthook import ExtDeprecationWarning
#warnings.simplefilter('ignore',ExtDeprecationWarning)
warnings.simplefilter('ignore',FutureWarning)
warnings.simplefilter('ignore',RuntimeWarning)

# Make sure that the camoco install is always frist
import sys
import os
import matplotlib
# Force matplotlib to not use any Xwindows backend.
matplotlib.use('Agg')

import camoco as co
import pandas as pd

import sys
import os
import glob
import argparse

from camoco.Tools import log 

pd.set_option('display.width',300)

# House Keeping commands
from camoco.cli.commands.build import build_cob,build_refgen,build_gont,build_GWAS
from camoco.cli.commands.remove import remove 
from camoco.cli.commands.list import list_command
# Functional and Validation commands
from camoco.cli.commands.health import cob_health
from camoco.cli.commands.simulateGWAS import simulateGWAS
from camoco.cli.commands.plotGWAS import plot_gwas
# Information Commands
from camoco.cli.commands.crossRef import crossref
from camoco.cli.commands.cistrans import cistrans
from camoco.cli.commands.snp2gene import snp2gene
from camoco.cli.commands.geneneighbors import geneneighbors
# Analysis Commands
from camoco.Overlap import Overlap

from camoco import cf

from argparse import ArgumentParser

if __name__ == '__main__':
    '''--------------------------
        Main Arguments 
    --------------------------'''
    parser = argparse.ArgumentParser(
        description=(
	"      ___           ___           ___           ___           ___           ___      \n"
	"     /  /\         /  /\         /__/\         /  /\         /  /\         /  /\     \n"
	"    /  /:/        /  /::\       |  |::\       /  /::\       /  /:/        /  /::\    \n"
	"   /  /:/        /  /:/\:\      |  |:|:\     /  /:/\:\     /  /:/        /  /:/\:\   \n"
	"  /  /:/  ___   /  /:/~/::\   __|__|:|\:\   /  /:/  \:\   /  /:/  ___   /  /:/  \:\  \n"
	" /__/:/  /  /\ /__/:/ /:/\:\ /__/::::| \:\ /__/:/ \__\:\ /__/:/  /  /\ /__/:/ \__\:\ \n"
	" \  \:\ /  /:/ \  \:\/:/__\/ \  \:\~~\__\/ \  \:\ /  /:/ \  \:\ /  /:/ \  \:\ /  /:/ \n"
	"  \  \:\  /:/   \  \::/       \  \:\        \  \:\  /:/   \  \:\  /:/   \  \:\  /:/  \n"
	"   \  \:\/:/     \  \:\        \  \:\        \  \:\/:/     \  \:\/:/     \  \:\/:/   \n"
	"    \  \::/       \  \:\        \  \:\        \  \::/       \  \::/       \  \::/    \n"
	"     \__\/         \__\/         \__\/         \__\/         \__\/         \__\/ \n\n\n" 
        "Camoco (Co-analysis of Molecular Components) inter-relates and co-analyzes different \n"
        "levels of genomic data. Namely it integrates genes present near and around GWAS loci\n"
        "using unbiased, functional information derived from co-expression networks."
	),
        formatter_class=argparse.RawDescriptionHelpFormatter,
        epilog= "\n".join([
            'version: {}'.format(co.__version__),
            'src:{}'.format(co.__file__),
            'Cache. Money. Corn.'
        ])
    )
    parser.add_argument(
        '--version',
        help='print the software version',
        default=False,
        action='store_true'
    )
    parser.add_argument(
        '--debug',
        action='store_true',
        default=False,
        help='Drop into ipdb when something bad happens.'
    )
    parser.add_argument(
        '--interactive',
        action='store_true',
        default=False,
        help='Initiate an ipdb session right before exiting.'
    )
    parser.add_argument(
        '--force',
        action='store_true',
        default=False,
        help='Overwrite output files from previous analyses.'
    )
    parser.set_defaults(func=lambda x: parser.print_help())
    subparsers = parser.add_subparsers(
        title='Camoco CLI program',
        metavar='Available Commands',
        description='Use --help with each command for more info',
    )


    # Also allow the help message to be printed using the help command
    helpcmd = subparsers.add_parser('help',help='Prints this help message')
    helpcmd.set_defaults(func=lambda x: parser.print_help())




    '''--------------------------
        Build GWAS
    --------------------------'''
    bldgwas = subparsers.add_parser(
        'build-gwas',
        help='build a GWAS dataset'
    )
    bldgwas.add_argument(
        'filename',
        help='Filename of csv/tsv.'
    )
    bldgwas.add_argument(
        'name',
        help='Name of GWAS dataset.'
    )
    bldgwas.add_argument(
        'description',
        help='Short description of dataset.'
    )
    bldgwas.add_argument(
        'refgen',
        help='Name of camoco RefGen dataset.'
    )
    bldgwas.add_argument(
        '--add-trait',
        default=None,
        help='Add a trait as a column to the CSV. This is '
        'useful if the file is only 1 trait and it is not '
        'present as a column.'
    )
    bldgwas.add_argument(
        '--trait-col',
        metavar='Trait',
        default='Trait',
        help='The name of the column in the table containing '
             'the trait name. '
    )
    bldgwas.add_argument(
        '--chrom-col',
        metavar='CHR',
        default='CHR',
        help=(
            'The name of the column containing the SNP chromosome (note: must '
            'correspond with chromosome names if RefGen)'
        )
    )
    bldgwas.add_argument(
        '--pos-col',
        metavar='POS',
        default='POS',
        help=(
            'The name of the column containing the SNP position (note: must '
            'correspond with the positions in RefGen)'
        )
    )
    bldgwas.add_argument(
        '--sep',
        metavar='\t',
        default='\t',
        help='Field Separator in CSV/TSV (default:\\t)'
    )
    bldgwas.add_argument(
        '--skip-traits',
        default=[],
        nargs='*',
        help=(
            'Skip these traits. Can provide as many traits as you like'
        )
    )
    bldgwas.add_argument(
        '--query',
        default=None,
        help=(
            'Before importing the traits from the dataframe, '
            'filter the columns via DataFrame.query() '
            'e.g. --query "pval < 0.05" would execute '
            'df.query("pval < 0.05") assuming the pval '
            'column exists. Useful for filtering data '
            'before importing the traits. '
            'See https://pandas.pydata.org/pandas-docs/stable/generated/pandas.DataFrame.query.html '
            'for more info.'
        )
    )
    bldgwas.set_defaults(func=build_GWAS)


    '''--------------------------
        Build GOnt
    --------------------------'''
    bldgo = subparsers.add_parser(
        'build-go',
        help='Build a Gene Ontology (GO)'
    )
    bldgo.add_argument(
        'filename',
        help='Gene-term map file'
    )
    bldgo.add_argument(
        'obo_filename',
        help='GO .obo filename'
    )
    bldgo.add_argument(
        'name',
        help='GO dataset name'
    )
    bldgo.add_argument(
        'description',
        help='short dataset description'
    )
    bldgo.add_argument(
        'refgen',
        help='Camoco reference Genome name.'
    )
    bldgo.add_argument(
        '--go-col',
        default=1,
        type=int,
        help='The GO Term ID column in Gene-term map file (default:1)'
    )
    bldgo.add_argument(
        '--id-col',
        default=0,
        type=int,
        help='The gene ID column in Gene-term map file (default:0)'
    )
    bldgo.add_argument(
        '--gene-term-header',
        default=False,
        action='store_true',
        help=(
            'Flag specifies if a header line is present in the gene-term map file'
            'If set, the build-command will skip a single (header) line in the file'
        )
    )
    bldgo.set_defaults(func=build_gont)
    

    '''--------------------------
        Build RefGen
    --------------------------'''
    bldref = subparsers.add_parser(
        'build-refgen',
        help='Build a Reference Genome.'
    )

    bldref.add_argument(
        'filename',
        help=(
            "The path to the GFF file."
        )
    )
    bldref.add_argument(
        'name',
        help=(
            "The name if the RefGen object to be stored in the core "
            "camoco database. "
        )
    )
    bldref.add_argument(
        'description',
        help=(
            "A short description of the RefGen for future reference "
        )
    )
    bldref.add_argument(
        'build',
        help=(   
            "A string designating the genome build, used for comparison"
            "operations, genes may share IDS but are different across build."
        )
    )
    bldref.add_argument(
        'organism',
        help=(
            'A short string describing the organims this RefGen is coming '
            'from. Again, is used in comparing equality among genes which '
            'may have the same id or name. '
        )
    )
    bldref.add_argument(
        '--chrom-feature',
        default='chromosome',
        metavar='chromosome',
        help=(
            'The name of the feature (in column 3) that designates a '
            'a chromosome. default: "chromosome"'
        )
    )
    bldref.add_argument(
        '--gene-feature',
        default='gene',
        metavar='gene',
        help=(
            'The name of the feature (in column 2) that designates a '
            'gene. These features will be the main object that the RefGen '
            'encompasses. default: "gene"'
        )
    )
    bldref.add_argument(
        '--ID-attr',
        default='ID',
        metavar='ID',
        help=(
            'The key in the attribute column which designates the ID or '
            'name of the feature. default: "ID"'
        )
    )
    bldref.add_argument(
        '--attr-split',
        default='=',
        metavar='=',
        help=(
            'The delimiter for keys and values in the attribute column. '
            'default: "="'
        )
    )
    bldref.set_defaults(func=build_refgen)

    '''--------------------------
        Build COB
    --------------------------'''
    bldcob = subparsers.add_parser(
        'build-cob',
        help='Build a Co-expression network.'
    )
    bldcob.add_argument(
        'filename',
        help="Path to CSV or TSV"
    )
    bldcob.add_argument(
        'name',
        help="Name of the network."
    )
    bldcob.add_argument(
        'description',
        help="Short description of network"
    )
    bldcob.add_argument(
        'refgen',
        type=str,
        help="Name of a pre-built RefGen. See build-refgen command"
    )
    bldcob.add_argument(
        '--rawtype',
        default='RNASEQ',
        help=(
            "Passing in a rawtype can help determine default parameters for "
            "standard data types such as MICROARRAY and RNASEQ. "
            "Specifies the fundamental datatype used to measure expression. "
            "During importation of the raw expression data, this value is "
            "used to make decisions in converting data to log-space. "
            "Options: (one of: 'RNASEQ' or 'MICROARRAY') "
        )
    )
    bldcob.add_argument(
        '--sep',
        default='\t',
        help="Field separators in the CSV or TSV, default=\\t"
    )
    bldcob.add_argument(
        '--index-col',
        type=str,
        default=0,
        help=(
            'If not None, this column will be set as the gene index '
            'column. Useful if there is a column name in the text file '
            'for gene names.'
        )
    )  
    bldcob.add_argument(
        '--max-val',
        type=float,
        help=(
            "This value is used to determine if any columns of the "
            "dataset have already been normalized. If any 'normailzed' "
            "values in an Accession column is larger than max_val, an "
            "exception is thown. max_val is determined by Expr.raw_type "
            "(default 100 for MicroArray and 1100 for RNASeq) but a "
            "max_val can be passed in to override these defaults. "
        )
    )  
    bldcob.add_argument(
        '--skip-quantile',
        action='store_true',
        default=False,
        help=(
            'Flag specifies whether or not to perform quantile normalization '
            'on expression data.'
        )
    ) 
    bldcob.add_argument(
        '--skip-quality-control',
        default=False,
        action='store_true', 
        help=(
            "A Flag indicating to skip quality control procedure "
            "on expression data. "
            "Default: False"
        )
    )  
    bldcob.add_argument(
        '--skip-normalization',
        default=False,
        action='store_true',
        help=(
            "Flag indicating that expression normaliztion should be skipped. "
            "Default: False"
        )
    )
    bldcob.add_argument(
        '--min-expr',
        default=0.001,
        type=float,
        help=(
            "Expression values (e.g. FPKM) under this threshold will be set "
            "to NaN and not used during correlation calculations. "
            "default: 0.001"
        )
    )  
    bldcob.add_argument(
        '--max-gene-missing-data',
        default=0.3,
        type=float,
        help=(   
            "Maximum percentage missing data a gene can have. Genes not "
            "meeting this criteria are removed from dataset."
            "default: 0.3"
        )
    )  
    bldcob.add_argument(
        '--max-accession-missing-data',
        default=0.5,
        type=float,
        help=(
            "maximum percentage missing data an accession (experiment) can"
            "have before it is removed."
            "Default: 0.5"
        )
    )  
    bldcob.add_argument(
        '--min-single-sample-expr',
        default=0.08,
        type=float,
        help=(
            "Genes that do not have a single accession having an expression "
            "value above this threshold are removed from analysis. These are "
            "likely presence/absence and will not have a strong coexpression "
            "pattern. Default: 0.08"
        )
    )  
    bldcob.add_argument(
        '--allow-non-membership',
        action='store_true',
        help=(
            "Flag indicating that feature (genes, metabolites, etc) should "
            "not be filtered out of the expression matrix because they are "
            "not present in the reference genome. This is useful for "
            "features, such as metabolites, that cannot be anchored to "
            "a RefGen. If true, features will be added to the RefGen with "
            "unknown coordinates"
        )
    )
    bldcob.add_argument(
        '--zscore_cutoff',
        default=3,
        metavar='(default: 3.0)',
        type=float,
        help=(
            "The zscore threshold used for edges in the network. Edges with "
            "z-scores under this value will not be used for thresholded "
            "calculations such as locality. Un-thresholded calculations, "
            "such as density will not be affected by this cutoff."
        )
    )
    bldcob.add_argument(
        '--dry-run',
        action='store_true',
        help='Dry run will only process the first 5000 genes'
    )
    bldcob.add_argument(
        '--non-interactive',
        action='store_true',
        help='If true, all the interactive steps/checks will be answered yes'
    )

    bldcob.set_defaults(func=build_cob)

    '''--------------------------
        list 
    --------------------------'''
    lstcmd = subparsers.add_parser(
        'list', 
        aliases=['ls'],
        help='List camoco datasets.' 
    )
    lstcmd.add_argument(
        'type', 
        nargs='?',
        default=None,
        help=(
            'The type of dataset you want to list. e.g. GWAS. '
            'Also accepts SQLite wildcards.'    
        )
    )
    lstcmd.add_argument(
        'name', 
        nargs='?',
        default=None,
        help = (
            'The name of the dataset you want to list. '
            'Also accepts SQLite wildcards.' 
        )
    )
    lstcmd.add_argument(
        '--names',
        action='store_true',
        help=(
            'If True, will return a list of dataset names'    
        )
    )
    lstcmd.add_argument(
        '--terms',
        action='store_true',
        help=(
            'Flag indicating to list the terms within the Ontology or GWAS' 
        )
    )
    lstcmd.set_defaults(func=list_command)

    ''' --------------------------
        rm 
    --------------------------'''
    removeCLI = subparsers.add_parser(
        'rm',
        help='Remove camoco dataset.'
    )
    removeCLI.add_argument('type',default=None)
    removeCLI.add_argument('name',default=None)
    removeCLI.add_argument('-f','--force',action='store_true',default=False)
    removeCLI.set_defaults(func=remove)

   #''' --------------------------
   #    plotGWAS 
   #--------------------------'''
   #plotgwas = subparsers.add_parser(
   #    'plotGWAS',
   #    help='Plot a GWAS dataset'
   #)
   ## Data Set Arguments
   #plotgwas.add_argument(
   #    'cob',
   #    help='The camoco network to use.'
   #)
   #plotgwas.add_argument(
   #    'gwas', 
   #    help='The camoco GWAS to use.'
   #)
   #plotgwas.add_argument(
   #    '--terms',
   #    nargs='*', 
   #    default=['all'],
   #    help='The term within the GWAS ontology to use. default: all'
   #)
   #plotgwas.add_argument(
   #    '--candidate-window-size',
   #    default=50000,
   #    type=int,
   #    help=(
   #        'The window size (in bp) for mapping effective SNPs to genes. '
   #        'default: 50000'
   #    )
   #)
   #plotgwas.add_argument(
   #    '--candidate-flank-limit',
   #    type=int,
   #    default=2,
   #    help=(
   #       'The number of flanking genes included in SNP to gene mapping. '
   #       'default: 2' 
   #    )
   #)
   #plotgwas.add_argument(
   #    '--snp2gene', 
   #    type=str,
   #    default='strongest',
   #    metavar='strongest',
   #    help=(
   #        'The SNP to gene mapping to use. Specifying multiple mappings '
   #        ' will result in multiple results. Must be one of: '
   #        "['effective','strongest'] "
   #        "default:'strongest' "
   #    )
   #)
   #plotgwas.add_argument(
   #    '--strongest-attr', 
   #    default='pval', 
   #    type=str,
   #    help=(
   #        'The locus attr used to determine which locus is the'
   #        'strongest locus. (defualt=pval).'
   #    )
   #)
   #plotgwas.add_argument(
   #    '--strongest-higher',
   #    action='store_false',
   #    default=True,
   #    help=(
   #        'Flag indicating the value in --strongest-attr is '
   #        'stronger if higher. Default behavior is to treat '
   #        'lower values as stronger (i.e. p-vals)'
   #    )
   #)
   #plotgwas.add_argument(
   #    '--out', 
   #    default=sys.stdout,
   #    help='OutputFile Name (default: Standard Out)'
   #)
   #plotgwas.set_defaults(func=plot_gwas) 

    ''' --------------------------
        Overlap (density/locality) 
    --------------------------'''
    overlapCLI = subparsers.add_parser(
        'overlap',
        help='Calculate network overlap among GWAS results. See --method for details.'
    )
    # Data Set Arguments
    overlapCLI.add_argument(
        'cob',
        help='The camoco network to use.'
    )
    overlapCLI.add_argument(
        'method',
        default='density',
        choices=['density','locality'],
        help='The metric used to evaulate overlap between Loci and the Network'
    )
    # keyword args
    overlapCLI.add_argument(
        '--genes',
        nargs='*',
        default=[None],
        help=(
            'Provide a [comma, space, semicolon] separated list of input genes.'    
        )
    )
    # Or get genes from terms
    overlapCLI.add_argument(
        '--gwas', 
        default=None,
        help='Calculate overlap on a Camoco GWAS dataset.'
    )

    overlapCLI.add_argument(
        '--go',
        default=None,
        help='Caluclate overlap among GO terms'
    )
    # supplementary args
    overlapCLI.add_argument(
        '--terms',
        nargs='*', 
        help=(
            'The term within the ontology to use. If all, '
            'terms in gwas will be iteratively analyzed. '
            '(default: all)'
        ),
        default=['all']
    )
    overlapCLI.add_argument(
        '--skip-terms',
        nargs='*',
        default=[],
        help=(
            'Terms specified here will be skipped.'    
        )
    )
    overlapCLI.add_argument(
        '--min-term-size',
        default=2,
        type=int,
        metavar='2',
        help='The minimum number of gene in a term to calculate overlap.'
    )
    overlapCLI.add_argument(
        '--max-term-size',
        default=None,
        type=int,
        metavar='None',
        help='The maximum number of genes in a term to calculate overlap.'
    )
    # SNP to gene mapping
    overlapCLI.add_argument(
        '--snp2gene', 
        type=str, 
        default='effective',
        choices=['strongest','effective'],
        metavar='effective',
        help=(
            'The SNP to gene mapping to use. If *effective*, SNPs with '
            'overlappingn windows will be collapsed into genomic intervals '
            'resulting in all genes within the intervals to be used. If '
            '*strongest* is specified, SNPs with lower values designated from '
            '--strongest-attr (e.g. pvals) will be dropped when SNPs have overlapping '
            'windows. Value must be in ["effective","strongest"]. ' 
            '(default: effective)'
        )
    )
    overlapCLI.add_argument(
        '--strongest-attr', 
        default='pval', 
        type=str, 
        metavar='pval',
        help=(
            'The locus attr used to determine which locus is the'
            'strongest locus. (defualt=pval).'
        )
    )
    overlapCLI.add_argument(
        '--strongest-higher',
        action='store_false',
        default=True,
        help=(
            'Flag indicating the value in --strongest-attr is'
            'stronger if higher. Default behavior is to treat'
            'lower values as stronger (i.e. p-vals)'
        )
    )
    overlapCLI.add_argument(
        '--candidate-window-size',
        type=int,
        default=1, 
        metavar=1,
        help=(
            'The window size, in bp, for mapping effective SNPs to genes. '
            '(default: 1)'
        )
    )
    overlapCLI.add_argument(
        '--candidate-flank-limit',
        type=int,
        default=0, 
        metavar=0,
        help=(
            'The number of flanking genes included in SNP to gene mapping. '
            'on each side of the locus. (default: 1)' 
        )
    )
    # Statistical Significance
    overlapCLI.add_argument(
        '--num-bootstraps', 
        default='auto',
        metavar='auto',
        help=(
            'The number of bootstraps to perform in order '
            'to estimate a null distribution. If auto the algorithm '
            'will bootstrap *until* the term is not significant at n=1000 '
            '*or* 1000 bootstraps have been performed. (default: auto)'
        )
    )

    # Output options 
    overlapCLI.add_argument(
        '--out', 
        default=None,
        help='OutputFile Name (default: Standard Out)'
    )
    overlapCLI.add_argument(
        '--dry-run',
        action='store_true',
        help='Do not perform the expensive computations'
    )
    overlapCLI.set_defaults(func=Overlap.from_CLI)


    ''' --------------------------
        Health
    --------------------------'''
    health = subparsers.add_parser(
        'health',
        help='Generate network health statistics'
    )
    health.add_argument(
        'cob',
        default=None
    )
    health.add_argument(
        '--out',
        default=None,
        help='Output file prefix.'
    )
    health.add_argument(
        '--refgen',
        default=None,
        help='Global Refgen is necessary for some stats'        
    )
    health.add_argument(
        '--go',
        default=None,
        help='Perform Gene Ontology Statistics. (default: None)'
    )
    health.add_argument(
        '--two-tailed',
        action='store_true',
        default=False,
        help='Include negative density values for GO statistics (default: False)'
    )
    health.add_argument(
        '--min-term-size',
        default=10,
        metavar=10,
        type=int,
        help='Minimum GO term size to consider in health (defualt: 10)'
    )
    health.add_argument(
        '--max-term-size',
        metavar=300,
        default=300,
        type=int,
        help='Maximum GO term size to consider in health (defualt: 300)'
    )
    health.add_argument(
        '--max-terms',
        default=None,
        type=int,
        help='If specified, limits the number of total GO terms tested (i.e. Subsampling)'
    )
    health.add_argument(
        '--edge-zscore-cutoff',
        default=None,
        type=int,
        metavar=3.0,
        help='Significance zscore cutoff for edges'
    )
    health.add_argument(
        '--num-bootstraps',
        default=100,
        type=int,
        metavar=100,
        help='Number of Bootstraps for pvalue stats.'
    )
    health.set_defaults(func=cob_health)


   #'''----------------------------------
   #    Cross Reference
   #-------------------------------------'''
   #crossref = subparsers.add_parser(
   #    'crossref',
   #    help=(
   #        'cross reference a set of networks'
   #        ' based on edge strength and go term density'
   #    )
   #)
   #crossref.add_argument(
   #    'cobs',
   #    default=None,
   #    nargs='*',
   #    help='A list of COBs which to cross reference'
   #)
   #crossref.add_argument(
   #    '--out',
   #    default=None,
   #    help='Output file name.'
   #)
   #crossref.set_defaults(func=crossref)

   #'''----------------------------------
   #    Simulate GWAS
   #-------------------------------------'''
   #simulateGWASCLI = subparsers.add_parser(
   #    'simulateGWAS',
   #    formatter_class=argparse.RawDescriptionHelpFormatter,
   #    help='Evaluate co-expression in simulated GWAS datasets',
   #    description=(
   #        'Simulate a GWAS dataset by adding candidate windows and \n'
   #        'flanking limits to terms within an Ontology. The ideal \n'
   #        'GWAS experiment is simulated by defining SNPs to having a 1:1 \n'
   #        'mapping to gene in a GO term. i.e. SNPs are directly generated \n'
   #        'using starting position of genes in GO term \n'
   #        'Noise is intruduced to the experiment \n'
   #        'using various sources of noise outlined below.\n'
   #        'Optional arguments below are applied in the following order: \n\n'
   #        '+ GO terms between min and max term size are extracted \n'
   #        '+ `--percent-SNPs-missed` of "SNPs" in GO term are removed. (False Negative SNPs) \n'
   #        '+ `--percent-fcr` of "SNPs" are replaced with random "SNPs" (False Positive SNPs) \n'
   #        '+ `--candidate-flank-limit` and `--candidate-window-size` introduce noise in \n'
   #        '   SNP to gene mapping.'
   #    )
   #)
   #simulateGWASCLI.add_argument(
   #    'cob',
   #    help='The Coexpression network you want to use.'
   #)
   #simulateGWASCLI.add_argument(
   #    'GOnt',
   #    help='The name of the GO Ontology you want to analyze.'
   #)
   #simulateGWASCLI.add_argument(
   #    '--terms',
   #    type=str,
   #    metavar='all',
   #    default = ['all'],
   #    nargs='*',
   #    help=("Ontology term to simulate. If 'all', which is the default, "
   #        "all terms will be tested."    
   #    )
   #)
   #simulateGWASCLI.add_argument(
   #    '--method',
   #    default='density',
   #    choices=['density','locality'],
   #    help='The metric used to evaluate overlap between the loci and the network'
   #)
   #simulateGWASCLI.add_argument(
   #    '--max-term-size',
   #    type=int,
   #    default=100,
   #    help='Maximum number of genes allowed in a GO term'
   #)
   #simulateGWASCLI.add_argument(
   #    '--min-term-size',
   #    type=int,
   #    default=50,
   #    help='Minimum number of genes allowed in a GO term'
   #)
   #simulateGWASCLI.add_argument(
   #    '--candidate-window-size',
   #    type=int,
   #    default=1, 
   #    metavar=1,
   #    help=(
   #        'The window size, in bp, for mapping effective SNPs to genes. '
   #    )
   #)
   #simulateGWASCLI.add_argument(
   #    '--candidate-flank-limit',
   #    type=int,
   #    default=0, 
   #    metavar=0,
   #    help=(
   #        'The number of flanking genes included in SNP to gene mapping. '
   #        'on each side of the locus.' 
   #    )
   #)
   #simulateGWASCLI.add_argument(
   #    '--num-bootstraps', 
   #    default='auto',
   #    metavar=100,
   #    help=(
   #        'The number of bootstraps to perform in order '
   #        'to estimate a null distribution.'
   #    )
   #)
   #simulateGWASCLI.add_argument(                                                  
   #    '--percent-fcr',                                                        
   #    type=int,                                                             
   #    default=0,                                                           
   #    metavar=0,                                                            
   #    help=(                                                                  
   #        '(FCR) False Cadidate Rate: '
   #        'The percentage of the term loci to be replaced with '              
   #        'random loci, simulating noise or false positives '                 
   #        'within the term. (Default: 0)'                                                   
   #    )                                                                       
   #)  
   #simulateGWASCLI.add_argument(
   #    '--percent-mcr',
   #    metavar=0,
   #    type=int,
   #    default=0,
   #    help=(
   #        '(MCR) Missing Candidate Rate: '
   #        'Percentage of input SNPs that were missed by '
   #        'simulated GWAS. e.g. if GWAS term SNPs ideally '
   #        'tag 10 GO genes and --percent-mcr=50 half '
   #        'of the genes in the GO term will be removed resulting '
   #        'in an input network of 5 genes. (Default=0)'
   #    )
   #)
   #simulateGWASCLI.add_argument(
   #    '--out', 
   #    default=None,
   #    help='OutputFile Name'
   #)
   #simulateGWASCLI.set_defaults(func=simulateGWAS.from_CLI)

   #'''----------------------------------
   #    CisTrans
   #-------------------------------------'''
   #cis_trans = subparsers.add_parser(
   #    'cistrans',
   #    help=(
   #        'compare the density of interactions between cis interactions'
   #        'to trans interactions'
   #    )
   #)
   #cis_trans.add_argument(
   #    'cob',
   #    default=None
   #)
   #cis_trans.add_argument(
   #    '--cis-distance',
   #    default=50000,
   #    type=int,
   #    help='Distance used to determine if two genes are in cis. (default: 50 kb)'
   #)
   #cis_trans.add_argument(
   #    '--out',
   #    default=None,
   #    help='Output file prefix'
   #)
   #cis_trans.set_defaults(func=cistrans)

    '''----------------------------------
        snp2gene
    -------------------------------------'''
    snp2gene_cmd = subparsers.add_parser(
        'snp2gene',
        help=(
            'Generate candidate genes and accompanying '
            'information from GWAS SNPs'    
        )
    )
    snp2gene_cmd.add_argument(
        'refgen',
        default=None,
        help=(
            'The reference genome to use. Can also provide the name of a '
            'COB object to use that refgen object'
        )
    )
    snp2gene_cmd.add_argument(
        'gwas',
        help='The camoco GWAS to generate candidates from.'
    )
    snp2gene_cmd.add_argument(
        '--terms',
        nargs='*',
        type=str,
        metavar='all',
        help=(
            'The term within the ontology to use. If all, '
            'terms in gwas will be iteratively analyzed. '
            'default: all'
        ),
        default=['all']
    )
    snp2gene_cmd.add_argument(
        '--candidate-window-size',
        type=int,
        nargs='*',
        default=[50000], 
        metavar=50000,
        help=(
            'The window size, in bp, for mapping effective SNPs to genes. '
            '(default: 50000)'
        )
    )
    snp2gene_cmd.add_argument(
        '--candidate-flank-limit',
        type=int,
        nargs='*',
        default=[1], 
        metavar=1,
        help=(
            'The number of flanking genes included in SNP to gene mapping. '
            'on each side of the locus. (default: 1)' 
        )
    )
    snp2gene_cmd.add_argument(
        '--snp2gene', 
        type=str, 
        default='effective',
        choices=['effective','strongest'],
        metavar='effective',
        help=(
            'The SNP to gene mapping to use. Specifying multiple mappings '
            ' will result in multiple results. Must be in: '
            '["effective","strongest"] '
            '(default: effective)'
        )
    )
    snp2gene_cmd.add_argument(
        '--strongest-attr', 
        default='pval', 
        type=str, 
        metavar='pval',
        help=(
            'The locus attr used to determine which locus is the '
            'strongest locus. (defualt=pval).'
        )
    )
    snp2gene_cmd.add_argument(
        '--include-parent-attrs',
        default=[],
        metavar='attrs',
        nargs='*',
        help=(
            'Include the specified parent locus attributes from '
            'the effective SNP.'
        )
    )
    snp2gene_cmd.add_argument(
        '--strongest-higher',
        action='store_false',
        default=True,
        help=(
            'Flag indicating the value in --strongest-attr is '
            'stronger if higher. Default behavior is to treat '
            'lower values as stronger (i.e. p-vals) '
        )
    )
    snp2gene_cmd.add_argument(
        '--gene-info',
        default=[],
        metavar='file',
        nargs='*',
        help=(
            'Include additional candidate gene info output from '
            'other tables. Tables will be joined on as many columns as '
            'possible. You can specify multiple files.'
        )
    )
    snp2gene_cmd.add_argument(
        '--out',
        default=sys.stdout
    )
    snp2gene_cmd.set_defaults(func=snp2gene)


    '''----------------------------------
      gene neighbors
    -------------------------------------'''
    geneneighbors_cmd = subparsers.add_parser(
        'neighbors',
        help=(
            'Generate significant gene neighbors '
            'from largest to smallest Z-score'    
        )
    )
    geneneighbors_cmd.add_argument(
        'cob',
        default=None
    )
    geneneighbors_cmd.add_argument(
        '--numneighbors',
        type=int,
        default=10, 
        metavar=10,
        help=(
            'The number of neighboring genes to return. '
            'on each side of the locus. (default: 10)' 
        )
    )
    geneneighbors_cmd.add_argument(
        '--zscore',
        type=int,
        default=3, 
        metavar=3,
        help=(
            'Z-score significance cutoff (default: 3). '
        )
    )
    geneneighbors_cmd.add_argument(
      '--out',
      default=sys.stdout
    )
    geneneighbors_cmd.set_defaults(func=geneneighbors)


    ''' -------------------------------------------------------------------
    # DO IT 
    -----------------------------------------------------------------------'''
    args = parser.parse_args()
    # Add debug options
    def main(args):
        '''
        returns a system code
        '''
        if args.version is True:
            print(co.__version__)
            return 0
        if args.debug is True:
            from IPython.core import ultratb
            sys.excepthook = ultratb.FormattedTB(
                mode='Verbose', color_scheme='Linux', call_pdb=1
            )
        # Execute the command
        try:
            return_value = args.func(args)
        except AttributeError as e:
            if args.debug:
                raise e
            sys.exit(e)
        except MemoryError as e:
            if args.debug:
                raise e
            log.warn(
                "Your computer ran out of memory! This analysis needs to "
                "be done on a machine with more memory (16Gb+ recommended)"        
            )
            sys.exit(e)
        except Exception as e:
            log.warn(
                "An error occured.\nrun with --debug or "
                "--interactive for more information\n\n" + str(e)
            )
            if args.debug:
                raise e
            else:
                sys.exit(1)
        if args.interactive is True:
            from IPython.core import ultratb
            sys.excepthook = ultratb.FormattedTB(
                mode='Verbose', color_scheme='Linux', call_pdb=1
            )
            from camoco.Exceptions import CamocoInteractive
            raise CamocoInteractive()
        # Return 
        return return_value
    return_value = main(args)
    sys.exit(return_value)
