<jupyter_start><jupyter_text>Evaluate mock community classification accuracyThe purpose of this notebook is to evaluate taxonomic classification accuracy of mock communities using different classification methods. Prepare the environment-----------------------First we'll import various functions that we'll need for generating the report.<jupyter_code>%matplotlib inline
from os.path import join, exists, expandvars
import pandas as pd
from IPython.display import display, Markdown
import seaborn.colors.xkcd_rgb as colors
from tax_credit.plotting_functions import (pointplot_from_data_frame,
                                           boxplot_from_data_frame,
                                           heatmap_from_data_frame,
                                           per_level_kruskal_wallis,
                                           beta_diversity_pcoa,
                                           average_distance_boxplots,
                                           rank_optimized_method_performance_by_dataset)
from tax_credit.eval_framework import (evaluate_results,
                                       method_by_dataset_a1,
                                       parameter_comparisons,
                                       merge_expected_and_observed_tables,
                                       filter_df,
                                       evaluate_unifrac)<jupyter_output><empty_output><jupyter_text>Configure local environment-specific values-------------------------------------------**This is the only cell that you will need to edit to generate basic reports locally.** After editing this cell, you can run all cells in this notebook to generate your analysis report. This will take a few minutes to run, as results are computed at multiple taxonomic levels.Values in this cell will not need to be changed, with the exception of ``project_dir``, to generate the default results contained within tax-credit. To analyze results separately from the tax-credit precomputed results, other variables in this cell will need to be set.<jupyter_code>## project_dir should be the directory where you've downloaded (or cloned) the 
## tax-credit repository. 
project_dir = expandvars("../../")

## expected_results_dir contains expected composition data in the structure
## expected_results_dir/<dataset name>/<reference name>/expected/
expected_results_dir = join(project_dir, "data/precomputed-results/", "mock-community")

## mock_results_fp designates the files to which summary results are written.
## If this file exists, it can be read in to generate results plots, instead
## of computing new scores.
mock_results_fp = join(expected_results_dir, 'mock_results.tsv')

## results_dirs should contain the directory or directories where
## results can be found. By default, this is the same location as expected 
## results included with the project. If other results should be included, 
## absolute paths to those directories should be added to this list.
results_dirs = [expected_results_dir]

## directory containing mock community data, e.g., feature table without taxonomy
mock_dir = join(project_dir, "data", "mock-community")

## Minimum number of times an OTU must be observed for it to be included in analyses. Edit this
## to analyze the effect of the minimum count on taxonomic results.
min_count = 1

## Define the range of taxonomic levels over which to compute accuracy scores.
## The default given below will compute order (level 2) through species (level 6)
taxonomy_level_range = range(2,7)


# we can save plots in this directory
outdir = join(expandvars("../../"), 'plots')
dataset_ids = ['mock-1', 'mock-2', 'mock-12', 'mock-13', 'mock-14', 'mock-15', 'mock-16', 'mock-18', 'mock-19', 
               'mock-20', 'mock-21', 'mock-22', 'mock-23']
method_ids = ['mindivlp', 'mindivlp_thresh', 'rdp', 'sortmerna', 'uclust', 'blast', 'blast+', 'naive-bayes', 'naive-bayes-bespoke', 'vsearch', 'mindivlp', 'mindivlp_thresh']
ref_ids = ['gg_13_8_otus']#, 'unite_20.11.2016_clean_fullITS']<jupyter_output><empty_output><jupyter_text>Find mock community pre-computed tables, expected tables, and "query" tables----------------------------------------------------------------------------Next we'll use the paths defined above to find all of the tables that will be compared. These include the *pre-computed result* tables (i.e., the ones that the new methods will be compared to), the *expected result* tables (i.e., the tables containing the known composition of the mock microbial communities), and the *query result* tables (i.e., the tables generated with the new method(s) that we want to compare to the *pre-computed result* tables).**Note**: if you have added additional methods to add, set `append=True`. If you are attempting to recompute pre-computed results, set `force=True`.This cell will take a few minutes to run if new results are being added, so hold onto your hat. If you are attempting to re-compute everything, it may take an hour or so, so go take a nap.<jupyter_code>mock_results = evaluate_results(results_dirs, 
                                expected_results_dir, 
                                mock_results_fp, 
                                mock_dir,
                                taxonomy_level_range=range(2,7), 
                                min_count=min_count,
                                taxa_to_keep=None, 
                                md_key='taxonomy', 
                                subsample=False,
                                per_seq_precision=True,
                                exclude=['other'],
                                dataset_ids=dataset_ids,
                                reference_ids=ref_ids,
                                method_ids=method_ids)<jupyter_output><empty_output><jupyter_text>Restrict analyses to a set of datasets or references: e.g., exclude taxonomy assignments made for purpose of reference database comparisons. This can be performed as shown below — alternatively, specific reference databases, datasets, methods, or parameters can be chosen by setting dataset_ids, reference_ids, method_ids, and parameter_ids in the evaluate_results command above.<jupyter_code>#mock_results = filter_df(mock_results, column_name='Parameters', values=['10000_14_0.01_8'], exclude=False)
mock_results = mock_results.reset_index(drop=True)<jupyter_output><empty_output><jupyter_text>Compute Unifrac Distances and Plot--------------------------Compute unifrac distances using evaluate_unifrac, which queries the 99 OTU taxonomy file for OTU IDs and computes unifrac using EMDUnifrac.<jupyter_code>taxonomy_fp = join(project_dir, "data/ref_dbs/gg_13_8_otus/99_otu_taxonomy.txt")
tree_fp = join(project_dir, "data/ref_dbs/gg_13_8_otus/trees/99_otus_unannotated.tree")  # May vary

df = evaluate_unifrac(
    taxonomy_fp,
    tree_fp,
    results_dirs, 
    expected_results_dir, 
    mock_results_fp, 
    mock_dir,
    min_count=min_count,
    taxa_to_keep=None, 
    md_key='taxonomy', 
    subsample=False,
    dataset_ids=dataset_ids,
    reference_ids=ref_ids,
    method_ids=method_ids)

df.to_csv(join(project_dir, "data/precomputed-results/mock-community/unifrac_data.tsv"), sep = "\t", index = False)<jupyter_output><empty_output><jupyter_text>Or load in unifrac distances.<jupyter_code>df = pd.read_csv(join(project_dir, "data/precomputed-results/mock-community/unifrac_data.tsv"), sep = "\t")<jupyter_output><empty_output><jupyter_text>Remove unwanted parameters from MinDivLP results.<jupyter_code>indicestoremove = list(df[ (df["Method"] == "mindivlp") & (df["Parameters"] != "10000_14_0.01_8") ].index.values)
indicestoremove += list(df[ (df["Method"] == "mindivlp_thresh") & (df["Parameters"] != "10000_14_0.01_8") ].index.values)
indicestoremove.sort()
new_results = df.drop(indicestoremove)
new_results[new_results["Method"] == "mindivlp"]["Parameters"].unique()
old_results = df
df = new_results
df["Parameters"].unique()<jupyter_output><empty_output><jupyter_text>Plot the data.<jupyter_code>color_palette={
    'expected': 'black', 'rdp': colors['baby shit green'], 'sortmerna': colors['macaroni and cheese'],
    'uclust': 'coral', 'blast': 'indigo', 'blast+': colors['electric purple'], 'naive-bayes': 'dodgerblue',
    'naive-bayes-bespoke': 'blue', 'vsearch': 'firebrick', 'mindivlp': colors['electric green'], 'mindivlp_thresh': colors['black']
}

y_vars = ["Unweighted Unifrac", "Weighted Unifrac"]
point = pointplot_from_data_frame(df, "Dataset", ["weighted-unifrac", "unweighted-unifrac"], 
                                  group_by="Reference", color_by="Method",
                                  color_palette=color_palette)
# Resize for better viewing
point['weighted-unifrac'].fig.set_figheight(6)
point['weighted-unifrac'].fig.set_figwidth(14)
point['unweighted-unifrac'].fig.set_figheight(6)
point['unweighted-unifrac'].fig.set_figwidth(14)
point['unweighted-unifrac'].fig
point['weighted-unifrac'].fig
for k, v in point.items():
    v.savefig(join(outdir, '{0}-lineplots.png'.format(k)))<jupyter_output><empty_output><jupyter_text>Compute and summarize precision, recall, and F-measure for mock communities----------------------------------------------------------------------------------------In this evaluation, we compute and summarize precision, recall, and F-measure of each result (pre-computed and query) based on the known composition of the mock communities. We then summarize the results in two ways: first with boxplots, and second with a table of the top methods based on their F-measures. **Higher scores = better accuracy**As a first step, we will evaluate **average** method performance at each taxonomic level for each method within each reference dataset type.**Note that, as parameter configurations can cause results to vary widely, average results are not a good representation of the "best" results. See [here](Optimized-method-performance) for results using optimized parameters for each method.**First we will define our [color palette](http://matplotlib.org/examples/color/named_colors.html) and the variables we want to plot. Via seaborn, we can apply the [xkcd crowdsourced color names](https://xkcd.com/color/rgb.txt). If that still doesn't match your hue, use hex codes.<jupyter_code>color_palette={
    'expected': 'black', 'rdp': colors['baby shit green'], 'sortmerna': colors['macaroni and cheese'],
    'uclust': 'coral', 'blast': 'indigo', 'blast+': colors['electric purple'], 'naive-bayes': 'dodgerblue',
    'naive-bayes-bespoke': 'blue', 'vsearch': 'firebrick', 'mindivlp': colors['electric green'], 'mindivlp_thresh': colors['burnt orange']
}

y_vars = ["Taxon Accuracy Rate", "Taxon Detection Rate", "L1 Error"]
# Optional for removing unwanted methods.

indicestoremove = list(mock_results[ (mock_results["Method"] == "mindivlp") & (mock_results["Parameters"] != "10000_14_0.01_8") ].index.values)
indicestoremove += list(mock_results[ (mock_results["Method"] == "mindivlp_thresh") & (mock_results["Parameters"] != "10000_14_0.01_8") ].index.values)
indicestoremove.sort()
new_results = mock_results.drop(indicestoremove)
new_results[new_results["Method"] == "mindivlp"]["Parameters"].unique()
old_results = mock_results
mock_results = new_results
mock_results["Parameters"].unique()
point = pointplot_from_data_frame(mock_results, "Level", y_vars, 
                                  group_by="Reference", color_by="Method",
                                  color_palette=color_palette)
for k, v in point.items():
    v.savefig(join(outdir, 'mock-{0}-lineplots.png'.format(k)))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1218860443065253, 'fraction_numerical': 0.010257736402034307, 'mean_word_length': 3.4030360531309296, 'pattern_counts': {'":': 0, '<': 39, '<?xml version=': 0, '>': 39, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 22, 'profanity_word_list.json': 2, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '639110', 'license': 'MIT', 'n_tokens_mistral': 3380, 'n_tokens_neox': 3187, 'n_words': 969, 'path': 'PythonCode/experiments/benchmark_vs_others/tax-credit-data/ipynb/mock-community/evaluate-classification-accuracy.ipynb', 'repo_name': 'KoslickiLab/DiversityOptimization', 'size': '11601'}