<jupyter_start><jupyter_text>Returns top mimicker and reverser drugs from a set of up and down genesets<jupyter_code>import re
import math
import time
import requests
import pandas as pd
import json
import scipy.stats as st
from IPython.display import display, IFrame, Markdown
import seaborn as sns
import matplotlib.pyplot as plt
%%appyter markdown
# {{title}}
L1000FWD_URL = 'https://maayanlab.cloud/L1000FWD'
pattern = ".*:(?P<drug>BRD-\w*)(-\d*){,3}:.*"
clustergrammer_url = 'https://amp.pharm.mssm.edu/clustergrammer/matrix_upload/'
def clustergrammer(df, name, clustergrammer_url, fignum=1, label="Clustergrammer"):
    clustergram_df = df.rename(columns={i:"Signature: %s"%i for i in df.columns}, index={i:"Drug: %s"%i for i in df.index})
    clustergram_df.to_csv(name, sep="\t")
    response = ''
    for i in range(5):
        try:
            res = requests.post(clustergrammer_url, files={'file': open(name, 'rb')})
            if not res.ok:
                response = res.text
                time.sleep(1)
            else:
                clustergrammer_url = res.text.replace("http:","https:")   
                break
        except Exception as e:
            response = e
            time.sleep(2)
    else:
        if type(response) == Exception:
            raise response
        else:
            raise Exception(response)
    display(IFrame(clustergrammer_url, width="1000", height="1000"))
    display(Markdown("**Figure %d** %s [Go to url](%s)"%(fignum, label, clustergrammer_url)))
cmap = sns.cubehelix_palette(50, hue=0.05, rot=0, light=1, dark=0)
def heatmap(df, filename, width=15, height=15):
#     fig = plt.figure(figsize=(width,height))
    cg = sns.clustermap(df, cmap=cmap, figsize=(width, height))
    cg.ax_row_dendrogram.set_visible(False)
    cg.ax_col_dendrogram.set_visible(False)
    display(cg)
    cg.savefig(filename)<jupyter_output><empty_output><jupyter_text>L1000FWD Drug MetadataGet drug metadata from L1000FWD [1]<jupyter_code>drug_meta = pd.read_csv("https://amp.pharm.mssm.edu/l1000fwd/download/Drugs_metadata.csv", index_col=0)
drug_meta.head(10)<jupyter_output><empty_output><jupyter_text>**Table 1** Drug metadata taken from L1000FWD  Empirical background distributionFetch precomputed background signatures. These were generated using CREEDS [2] signatures<jupyter_code>empirical_mimickers_dist = pd.read_csv("https://appyters.maayanlab.cloud/storage/L1000FWDConsensus/creeds_drug_stats/CREEDS_mimickers_stat_new.tsv", sep="\t", index_col=0)
empirical_mimickers_dist.head(10)<jupyter_output><empty_output><jupyter_text>**Table 2** Empirical mean and standard deviation for mimicker drugs generated using CREEDS signatures<jupyter_code>empirical_reversers_dist = pd.read_csv("https://appyters.maayanlab.cloud/storage/L1000FWDConsensus/creeds_drug_stats/CREEDS_reversers_stat_new.tsv", sep="\t", index_col=0)
empirical_reversers_dist.head(10)<jupyter_output><empty_output><jupyter_text>**Table 3** Empirical mean and standard deviation for reverser drugs generated using CREEDS signatures  Read Input<jupyter_code>%%appyter code_exec

{% set up_gene_sets = FileField(
    name='up_gene_sets',
    label='up gene sets',
    default='up_gene_sets.gmt',
    section="PRIMARY",
    examples={
        'up_gene_sets.gmt': 'https://appyters.maayanlab.cloud/storage/L1000FWDConsensus/sample_data/up_diseases'
    }
) %}

{% set down_gene_sets = FileField(
    name='down_gene_sets',
    label='down gene sets',
    default='down_gene_sets.gmt',
    section="PRIMARY",
    examples={
        'down_gene_sets.gmt': 'https://appyters.maayanlab.cloud/storage/L1000FWDConsensus/sample_data/down_diseases'
    }
) %}

up_gene_sets = {{ up_gene_sets }}
down_gene_sets = {{ down_gene_sets }}
signatures = {}
with open(up_gene_sets) as upfile:
    for line in upfile:
        unpacked = line.strip().split("\t\t")
        if not len(unpacked) == 2:
            raise ValueError("GMT is not formatted properly, please consult the README of the appyter for proper formatting")
        sigid, geneset_str = unpacked
        geneset = geneset_str.split("\t")
        signatures[sigid] = {
            "up_genes": [i.split(",")[0] for i in geneset]
        }
with open(down_gene_sets) as downfile:
    for line in downfile:
        unpacked = line.strip().split("\t\t")
        if not len(unpacked) == 2:
            raise ValueError("GMT is not formatted properly, please consult the README of the appyter for proper formatting")
        sigid, geneset_str = unpacked
        geneset = geneset_str.split("\t")
        if sigid not in signatures:
            raise ValueError("%s did not match any of the up signatures, make sure that the signature names are the same for both up and down genes")
        else:
            signatures[sigid]["down_genes"] = [i.split(",")[0] for i in geneset]
input_sigs = pd.DataFrame.from_dict(signatures, orient="index")
input_sigs.head(10)<jupyter_output><empty_output><jupyter_text>**Table 4** Input signatures  L1000FWDPerform signature search using L1000FWD (https://amp.pharm.mssm.edu/l1000fwd/)<jupyter_code>num_sigs = len(signatures)
print("Total number of signatures: %d"% num_sigs)
sig_search_endpoint = "/sig_search"
count = 1
response = ''
for k,v in signatures.items():    
    print("Querying (%d/%d): %s"%(count, num_sigs, k), end="\r", flush=True)
    count += 1
    payload = v
    for tries in range(5):
        try:
            res = requests.post(L1000FWD_URL + sig_search_endpoint, json=payload)
            if res.ok:
                result_id = res.json()["result_id"]
                signatures[k]["result_id"] = result_id
                time.sleep(0.1)
                break
            else:
                response = res.text
                time.sleep(0.2)
        except Exception as e:
            response = e
            time.sleep(1)
    else:
        if type(response) == Exception:
            raise response
        else:
            raise Exception(response)
topn_endpoint = "/result/topn/"
count = 1
response = ''
for k,v in signatures.items():    
    print("Querying (%d/%d): %s"%(count, num_sigs, k), end="\r", flush=True)
    count += 1
    endpoint = topn_endpoint + v["result_id"]
    for tries in range(5):
        try:
            res = requests.get(L1000FWD_URL + endpoint)
            if res.ok:
                results = res.json()
                signatures[k].update(results)
                time.sleep(0.1)
                break
            else:
                response = res.text
                time.sleep(0.2)
        except Exception as e:
            response = e
            time.sleep(1)
    else:
        if type(response) == Exception:
            raise response
        else:
            raise Exception(response)
with open("signatures.json", "w") as o:
    o.write(json.dumps(signatures, indent=4))<jupyter_output><empty_output><jupyter_text>[Download signature search results](./signatures.json)  User defined parameters<jupyter_code>%%appyter code_exec
alpha = {{FloatField(name='alpha', label='p-value cutoff', default=0.05, section='PRIMARY')}}
drug_perc = {{FloatField(name='drug_perc',
                         label='drug percentage',
                         description='Minimum percentage cutoff for drug instances in signatures',
                         default=0.1, section='PRIMARY')}}
width = {{FloatField(name='width', label='image width', default=10, section='PRIMARY')}}
height = {{FloatField(name='height', label='image height', default=10, section='PRIMARY')}}<jupyter_output><empty_output><jupyter_text>Mimickers<jupyter_code>mimickers_df = pd.DataFrame(index=signatures.keys())
for k,v in signatures.items():
    for sig in v["similar"]:
        sig_id = sig["sig_id"]
        m = re.match(pattern, sig_id)
        pert_id = m.group("drug")
        drug_name = drug_meta.loc[pert_id, "pert_iname"]
        drug = pert_id if pert_id == drug_name else "%s-%s" %(pert_id, drug_name)
        if drug not in mimickers_df.columns:
            mimickers_df[drug] = 0
        mimickers_df.at[k, drug] += 1
mimickers_df = mimickers_df.transpose()

# Get only drugs that appear on a percentage of the queries (defined by drug_perc)
mimickers_df = mimickers_df[(mimickers_df > 0).sum(1) > (num_sigs * drug_perc)]
mimickers_df.to_csv("mimickers.tsv", sep="\t")
mimickers_df.head(10)<jupyter_output><empty_output><jupyter_text>**Table 5** Mimicker drug counts ([download](./mimickers.tsv))<jupyter_code>mimickers_df.shape<jupyter_output><empty_output><jupyter_text>Consensus mimicker drugs<jupyter_code># Background signatures
empirical_mimickers_stat = empirical_mimickers_dist.loc[mimickers_df.index]
top_mimicker_drugs = mimickers_df.sum(1).to_frame(name="counts")
top_mimicker_drugs["number of signature queries"] = (mimickers_df > 0).sum(1)
top_mimicker_drugs["sample mean"] = mimickers_df.mean(1)
top_mimicker_drugs["sample stdev"] = mimickers_df.std(1)
top_mimicker_drugs["empirical mean (creeds)"] = empirical_mimickers_stat["mean"]
top_mimicker_drugs["empirical stdev (creeds)"] = empirical_mimickers_stat["stdev"]<jupyter_output><empty_output><jupyter_text>z-test$z=\frac{\bar{x}-\mu}{\frac{\sigma}{\sqrt{n}}}$Where:- $\bar{x}$ is the sample mean for a drug- $\mu$ is the computed empirical mean of the drugs using CREEDS signatures- $\sigma$ is the computed empirical standard deviation of the drugs using CREEDS signatures- $n$ is the number of samples<jupyter_code># Compute zstat and p value
top_mimicker_drugs["zstat"] = (top_mimicker_drugs["sample mean"] - empirical_mimickers_stat["mean"]) / (empirical_mimickers_stat["stdev"]/math.sqrt(num_sigs))
top_mimicker_drugs["p value"] = top_mimicker_drugs.zstat.apply(lambda x: 1-st.norm.cdf(x))

#Filter by p-value
top_mimicker_drugs = top_mimicker_drugs[top_mimicker_drugs["p value"]<alpha].sort_values(by=["p value"])

# Save to tsv
top_mimicker_drugs.to_csv("consensus_mimicker_drugs.tsv", sep="\t")
top_mimicker_drugs.head(10)<jupyter_output><empty_output><jupyter_text>**Table 6** Consensus mimicker drugs ([download](./consensus_mimicker_drugs.tsv)) Table shape<jupyter_code>top_mimicker_drugs.shape
top_mimickers_df = mimickers_df.loc[top_mimicker_drugs.index]
if (top_mimickers_df.shape[1] > 0):
    clustergrammer(top_mimickers_df, "consensus_mimickers_clust.tsv", clustergrammer_url, 1, "Consensus mimicker drugs clustergrammer")
else:
    print("No drugs found")
if (top_mimickers_df.shape[1] > 0):
    heatmap(top_mimickers_df, "consensus_mimickers.svg", width, height)
else:
    print("No drugs found")<jupyter_output><empty_output><jupyter_text>**Figure 2** Consensus mimickers heatmap [download figure](./consensus_mimickers.svg)  Reversers<jupyter_code>reversers_df = pd.DataFrame(index=signatures.keys())
for k,v in signatures.items():
    for sig in v["opposite"]:
        sig_id = sig["sig_id"]
        m = re.match(pattern, sig_id)
        pert_id = m.group("drug")
        drug_name = drug_meta.loc[pert_id, "pert_iname"]
        drug = pert_id if pert_id == drug_name else "%s-%s" %(pert_id, drug_name)
        if drug not in reversers_df.columns:
            reversers_df[drug] = 0
        reversers_df.at[k, drug] += 1
reversers_df = reversers_df.transpose()

# Get only drugs that appear on a percentage of the queries (defined by drug_perc)
reversers_df = reversers_df[(reversers_df > 0).sum(1) > (num_sigs* drug_perc)]
reversers_df.to_csv("reversers.tsv", sep="\t")
reversers_df.head(10)<jupyter_output><empty_output><jupyter_text>**Table 7** Reverser drug counts ([download](./reversers.tsv))  Consensus reverser drugs<jupyter_code># Get background distribution
empirical_reversers_stat = empirical_reversers_dist.loc[reversers_df.index]
top_reverser_drugs = reversers_df.sum(1).to_frame(name="counts")
top_reverser_drugs["number of signature queries"] = (reversers_df > 0).sum(1)
top_reverser_drugs["sample mean"] = reversers_df.mean(1)
top_reverser_drugs["sample stdev"] = reversers_df.std(1)
top_reverser_drugs["empirical mean (creeds)"] = empirical_reversers_stat["mean"]
top_reverser_drugs["empirical stdev (creeds)"] = empirical_reversers_stat["stdev"]<jupyter_output><empty_output><jupyter_text>z-test$z=\frac{\bar{x}-\mu}{\frac{\sigma}{\sqrt{n}}}$Where:- $\bar{x}$ is the sample signatures for a drug- $\mu$ is the computed empirical mean of the drugs using CREEDS signatures- $\sigma$ is the computed empirical standard deviation of the drugs using CREEDS signatures- $n$ is the number of signatures<jupyter_code># Compute zstat and p-value
top_reverser_drugs["zstat"] = (top_reverser_drugs["sample mean"] - empirical_reversers_stat["mean"]) / (empirical_reversers_stat["stdev"]/math.sqrt(num_sigs))
top_reverser_drugs["p value"] = top_reverser_drugs.zstat.apply(lambda x: 1-st.norm.cdf(x))

# Filter by p-value
top_reverser_drugs = top_reverser_drugs[top_reverser_drugs["p value"]<alpha].sort_values(by=["p value"])

# Save tsv
top_reverser_drugs.to_csv("consensus_reverser_drugs.tsv", sep="\t")
top_reverser_drugs.head(10)<jupyter_output><empty_output><jupyter_text>**Table 8** Consensus reverser drugs ([download](./consensus_reverser_drugs.tsv)) Table shape<jupyter_code>top_reverser_drugs.shape
top_reversers_df = reversers_df.loc[top_reverser_drugs.index]
if (top_reversers_df.shape[1] > 0):
    clustergrammer(top_reversers_df, "consensus_reversers_clust.tsv", clustergrammer_url, 3, "Consensus reverser drugs clustergrammer")
else:
    print("No drugs found")
if (top_reversers_df.shape[1] > 0):
    heatmap(top_reversers_df, "consensus_reversers.svg", width, height)
else:
    print("No drugs found")<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12743871015239638, 'fraction_numerical': 0.012883751748509167, 'mean_word_length': 4.267157813105855, 'pattern_counts': {'":': 1, '<': 68, '<?xml version=': 0, '>': 76, 'https://': 8, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '16', 'dir': 'jupyter-structured-clean-dedup', 'id': '63731', 'license': 'CC0-1.0', 'n_tokens_mistral': 4938, 'n_tokens_neox': 4674, 'n_words': 934, 'path': 'appyters/L1000FWD_Consensus_Drugs/L1000FWD Consensus Drugs.ipynb', 'repo_name': 'shui02/appyter-catalog', 'size': '13583'}