<jupyter_start><jupyter_text>Cell Line AnalysisWe sought to validate the Ras classifier trained on TCGA pan-cancer data by generating predictions on cell line data. A good classifier should generalize to predicting Ras status in other samples. We apply the classifier on two datasets:1. [GSE94937](https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE94937) from [Kim et al. 2017](http://doi.org/10.1016/j.cels.2017.08.002)  * Illumina NextSeq 5000 of Human Small Airway Epithelial Cells expressing KRAS G12V and wild-type KRAS2. [Cancer Cell Line Encyclopedia (CCLE)](https://software.broadinstitute.org/software/cprg/?q=node/11) Gene Expression data.  * 737 cell lines with matching gene expression and mutation calls  * Pharmacologic profiling of 24 drugs over 504 cell lines  These data were accessed via publicly available resources with help from links in the [UCSD-CCAL Onco-GPS github repository](https://github.com/UCSD-CCAL/onco-gps-paper-analysis)<jupyter_code>import os
import numpy as np
import pandas as pd
from decimal import Decimal
from scipy.stats import ttest_ind
from statsmodels.stats.proportion import proportions_chisquare
from sklearn.preprocessing import StandardScaler
from Bio.SeqUtils import IUPACData

import matplotlib.pyplot as plt
import seaborn as sns
import plotnine as gg
# Store protein change dictionary
aa = IUPACData.protein_letters_1to3_extended
%matplotlib inline<jupyter_output><empty_output><jupyter_text>Load PI3K_gain Classifier Coefficients<jupyter_code>classifier_file = os.path.join('..', 'classifiers', 'ERBB2_PIK3CA_KRAS_AKT1', 'classifier_summary.txt')
with open(classifier_file) as class_fh:
    for line in class_fh:
        line = line.strip().split('\t')
        if line[0] == 'Coefficients:':
            all_coef_df = pd.read_table(os.path.join('..', line[1]), index_col=0)

# Only non-zero coefficients contribute to model performance
coef_df = all_coef_df[all_coef_df['abs'] > 0]
print(coef_df.shape)
coef_df.head(10)<jupyter_output>(301, 3)<jupyter_text>Part 2: CCLENote - This data was also retrieved from the Onco-GPS paper analysis repository<jupyter_code>#ccle_file_name = os.path.join('..', '..', 'onco-gps-paper-analysis', 'data',
#                              'rpkm__gene_x_ccle_cellline.gct')
ccle_file_name = os.path.join('..','data','ccle_rnaseq_genes_rpkm_20180929_mod.gct')
ccle_df = pd.read_csv(ccle_file_name, skiprows=2, index_col=0)
print(ccle_df.shape)
ccle_df = ccle_df.drop_duplicates(subset='Description',keep = 'first')
ccle_df.head()
# Subset to common genes in the classifier and CCLE data
common_genes = list(set(coef_df['feature']) & set(ccle_df.index))
common_ccle_coef = coef_df[coef_df['feature'].isin(common_genes)]
print(common_ccle_coef.shape)
common_ccle_coef.head()
ccle_df = ccle_df.loc[common_ccle_coef['feature'], ccle_df.columns[1:]]
ccle_df.head()
scaled_fit = StandardScaler().fit(ccle_df.T)
ccle_df = pd.DataFrame(scaled_fit.transform(ccle_df.T),
                            index=ccle_df.columns,
                            columns=ccle_df.index)

print(ccle_df.shape)
ccle_df.head()
ccle_df = ccle_df.T
print(ccle_df.shape)
#ccle_df = ccle_df2.drop_duplicates(subset = ccle_df2.index, keep = 'first')
ccle_df.head()
# Get the weights ready for applying the classifier
apply_weights = pd.DataFrame(common_ccle_coef['weight'])
apply_weights.index = common_ccle_coef.feature
print(apply_weights.shape)
apply_weights.head()
# Apply a logit transform [y = 1/(1+e^(-wX))] to output probabilities
result_ccle = apply_weights.T.dot(ccle_df)
result_ccle = 1 / (1 + np.exp(-1 * result_ccle))
# Distribution of predictions of the Ras Classifier applied to CCLE data
result_ccle.T.hist();
# Load CCLE Mutation Data
#ccle_mut_file_name = os.path.join('..', '..', 'onco-gps-paper-analysis', 'data', 
#                                  'mutation__gene_x_ccle_cellline.gct')
ccle_mut_file_name = os.path.join('..','data','CCLE_MUT_CNA_AMP_DEL_binary_Revealer.gct')
ccle_all_mut_df = pd.read_table(ccle_mut_file_name, skiprows=2, index_col=0)
print(ccle_all_mut_df.shape)
ccle_all_mut_df.head()
# Load CCLE Variant Data
#ccle_maf_file = 'https://data.broadinstitute.org/ccle/CCLE_DepMap_18Q1_maf_20180207.txt'
ccle_maf_file = os.path.join('..','data','CCLE_DepMap_18Q1_maf_20180207.txt')
ccle_maf_df = pd.read_table(ccle_maf_file, index_col=15)
print(ccle_maf_df.shape)
ccle_maf_df.head(2)
# Identify all cell lines with mutations in pi3k_gain genes, also subset pten mutant samples
pi3k_genes = ['AKT1_MUT', 'PIK3CA_MUT', 'KRAS_MUT', 'ERBB2_MUT']
pi3k_status = ccle_all_mut_df.loc[pi3k_genes, :].T.apply(max, axis=1)

# BRAF mutations do not contribute to Ras status in this case
ccle_mut_df = (
    ccle_all_mut_df.loc[pi3k_genes + ['BRAF_MUT'], :].T
    .assign(pi3k_status=pi3k_status).drop(['Description'])
    )
print(ccle_mut_df.shape)
ccle_mut_df.head(5)
# Join classifier scores with mutation status
ccle_full_df = ccle_mut_df.join(result_ccle.T).dropna()
ccle_full_df = ccle_full_df.assign(sample_name = ccle_full_df.index)
ccle_full_df = ccle_full_df.sort_values(by='weight', ascending=False)
ccle_full_df.index.name = 'cell_line'
print(ccle_full_df.shape)
ccle_full_df.head(5)
# Write CCLE Scores to file
ccle_scores_file = os.path.join('..', 'results', 'ccle_pi3k_gain_BRAF_classifier_scores.tsv')
ccle_full_df.to_csv(ccle_scores_file, sep='\t')
# Use Seaborn for the 2nd plot
sns.set(style="whitegrid")
sns.set_context("paper", rc={"font.size":11, "axes.titlesize":11, "axes.labelsize":16,
                             'xtick.labelsize':11, 'ytick.labelsize':11})<jupyter_output><empty_output><jupyter_text>Perform a t-test on classifier weights across groups<jupyter_code># pi3k mutant vs. pi3k wildtype
pi3k_mutant = ccle_full_df[ccle_full_df['pi3k_status'] == 1]
pi3k_wt = ccle_full_df[ccle_full_df['pi3k_status'] == 0]

# Also interested in BRAF status within Ras wildtype samples
BRAF_mutant = pi3k_wt[pi3k_wt['BRAF_MUT'] == 1]
BRAF_wt = pi3k_wt[pi3k_wt['BRAF_MUT'] == 0]

# Also interested in PTEN status within PI3k gain_mutant samples
BRAF_mutant_2 = pi3k_mutant[pi3k_mutant['BRAF_MUT'] == 1]
BRAF_wt_2 = pi3k_mutant[pi3k_mutant['BRAF_MUT'] == 0]

# Output t-test results
t_results_pi3k = ttest_ind(a = pi3k_mutant['weight'],
                          b = pi3k_wt['weight'], equal_var = False)
print('pi3k Status:')
print(t_results_pi3k)

t_results_BRAF = ttest_ind(a = BRAF_mutant['weight'],
                           b = BRAF_wt['weight'], equal_var = False)
print('\nBRAF Status in PI3K Wild-Type Samples:')
print(t_results_BRAF)

t_results_BRAF_2 = ttest_ind(a = BRAF_mutant_2['weight'],
                           b = BRAF_wt_2['weight'], equal_var = False)
print('\nBRAF Status in PI3K mutant-Type Samples:')
print(t_results_BRAF_2)

# Plot Results
x1, x2 = 0, 1
x3, x4 = -0.2, 0.2
x5, x6 = 0.8, 1.2
y1, y2,y3, h = 1.17, 1.0,1.02, 0.03


plt.rcParams['figure.figsize']=(3.5, 4)
ax = sns.boxplot(x="pi3k_status", y="weight", data=ccle_full_df,
                 hue='BRAF_MUT', palette = {0: "whitesmoke", 1: 'gainsboro'},
                 fliersize=0)
ax = sns.stripplot(x='pi3k_status', y='weight', hue='BRAF_MUT',
                   data=ccle_full_df, 
                   dodge=True, edgecolor='gray',
                   palette = {1: "seagreen", 0: 'goldenrod'},
                   jitter=0.25, size=2, alpha=0.65)
handles, labels = ax.get_legend_handles_labels()
l = plt.legend(handles[2:4], ['Wild-Type', 'Mutant'], bbox_to_anchor=(.63, 0.2), loc=2, borderaxespad=0.)
l.set_title("BRAF")
ax.axes.set_ylim(0, 1.3)
ax.set_yticklabels([0, 0.2, 0.4, 0.6, 0.8, 1, ''])
ax.set_xticklabels(['Pi3k Wild-Type', 'Pi3k Mutant'])
ax.set_ylabel('Pi3k Classifier Score')
ax.set_xlabel('CCLE Data')
ax.legend
plt.axhline(0.5, color='black', linestyle='dashed', linewidth=1)


# Add PI3K T-Test Results
plt.plot([x1, x1, x2, x2], [y1, y1+h, y1+h, y1], lw=1.2, c='black')
plt.text(.5, y1+h, "{:.2E}".format(Decimal(t_results_pi3k.pvalue)),
         ha='center', va='bottom', color="black")

# Add BRAF t-test results
plt.plot([x3, x3, x4, x4], [y2, y2+h, y2+h, y2], lw=1.2, c='black')
plt.text(0, y2+h, "{:.2E}".format(Decimal(t_results_BRAF.pvalue)),
         ha='center', va='bottom', color="black")

# Add BRAF t-test results
plt.plot([x5, x5, x6, x6], [y3, y3+h, y3+h, y3], lw=1.2, c='black')
plt.text(1, y3+h, "{:.2E}".format(Decimal(t_results_BRAF_2.pvalue)),
         ha='center', va='bottom', color="black")

plt.tight_layout()
ccle_fig_file = os.path.join('..', 'figures', 'cell_line', 'ccle_pi3k_wild_mut_BRAF_predictions.pdf')
plt.savefig(ccle_fig_file)<jupyter_output><empty_output><jupyter_text>What percentage of correct classifications in CCLE data?<jupyter_code># Assign a label to what the predictions are given classifier scores
ccle_full_df = ccle_full_df.assign(predictions = 'wild-type')
ccle_full_df.loc[ccle_full_df['weight'] > 0.5, 'predictions'] = 'mutant'
# Stratify cell lines based on predictions and ground truth status
positive_pi3k_predictions_ccle = ccle_full_df[ccle_full_df['weight'] > 0.5]
negative_pi3k_predictions_ccle = ccle_full_df[ccle_full_df['weight'] <= 0.5]

positive_pi3k_lines_ccle = ccle_full_df[ccle_full_df['pi3k_status'] == 1]
negative_pi3k_lines_ccle = ccle_full_df[ccle_full_df['pi3k_status'] == 0]
# Of wild-type Ras cell lines, how many are predicted correctly?
# True Negative Rate, Specificity
negative_pi3k_lines_ccle['predictions'].value_counts()
# Of mutated Ras cell lines, how many are predicted correctly?
# True Positive Rate (TPR), Recall, Sensitivity
positive_pi3k_lines_ccle['predictions'].value_counts()
# Of the wild-type predictions, how many are actually wild-type?
# Negative Predictive Value (NPV)
neg_ccle_results = negative_pi3k_predictions_ccle['pi3k_status'].value_counts()
true_neg = neg_ccle_results[0]
predicted_condition_neg = neg_ccle_results.sum()

print('{} out of {} Pi3k wild-type predictions '
      'are true ({:.1f}%)'.format(true_neg, predicted_condition_neg,
                                  true_neg * 100 / predicted_condition_neg))
# Of the mutated predictions, how many are actually mutated?
# Positive Predictive Value (PPV) -or- precision
pos_ccle_results = positive_pi3k_predictions_ccle['pi3k_status'].value_counts()
false_pos, true_pos = pos_ccle_results
predicted_condition_pos = pos_ccle_results.sum()

print('{} out of {} Pi3k mutation predictions '
      'are true ({:.1f}%)'.format(true_pos, predicted_condition_pos,
                                  true_pos * 100 / predicted_condition_pos))
total_correct = true_pos + true_neg
print('{} of {} Total cell lines '
      'predicted correctly ({:.1f}%)'.format(total_correct, ccle_full_df.shape[0],
                                             total_correct * 100 / ccle_full_df.shape[0]))
# Of the False positives, how many are BRAF mutant?
wt_pi3k_braf_ccle = positive_pi3k_predictions_ccle[positive_pi3k_predictions_ccle['pi3k_status'] == 0]
braf_neg, braf_pos = wt_pi3k_braf_ccle['BRAF_MUT'].value_counts()

print('{} of {} total false positives '
      'have BRAF mutations ({:.1f}%)'.format(braf_pos, false_pos,
                                             braf_pos * 100 / false_pos))
# If include BRAF mutations, how many correct
correct_braf = wt_pi3k_braf_ccle['BRAF_MUT'].value_counts()[1]
true_pos_with_braf = true_pos + correct_braf
print('Including BRAF mutants, {} of {} PI3K mutation predictions '
      'have pi3k pathway mutations ({:.1f}%)'.format(true_pos_with_braf,
                                                    predicted_condition_pos,
                                                    true_pos_with_braf * 100 / predicted_condition_pos))
print('Of the false positives, there are {} BRAF mutated cell lines '
      'and {} BRAF wild-type cell lines'.format(braf_pos, braf_neg))
total_braf_wildtype, total_braf_mut = ccle_full_df['BRAF_MUT'].value_counts()
print('In all of CCLE, there are {} BRAF mutated cell lines '
      'and {} BRAF wild-type cell lines'.format(total_braf_mut, total_braf_wildtype))<jupyter_output>In all of CCLE, there are 87 BRAF mutated cell lines and 826 BRAF wild-type cell lines<jupyter_text>What Ras mutations are identified in False Negatives?**Variant level predictions in the CCLE**<jupyter_code># How many MAF tumors also have CCLE and mutation data?
len(set(ccle_full_df.index).intersection(set(ccle_maf_df.index)))
# Of the false negaves, what RAS mutations do they harbor?
false_negatives_df = negative_pi3k_predictions_ccle.query('pi3k_status == 1')
common_neg_ccle_samples = set(ccle_maf_df.index).intersection(set(false_negatives_df.index))
false_neg_maf_df = ccle_maf_df.loc[list(common_neg_ccle_samples), :]
# What about true positives?
true_positives_df = positive_pi3k_predictions_ccle[positive_pi3k_predictions_ccle['pi3k_status'] == 1]
common_pos_ccle_samples = set(ccle_maf_df.index).intersection(set(true_positives_df.index))
true_pos_maf = ccle_maf_df.loc[list(common_pos_ccle_samples), :]
# Subset to only Ras genes 'AKT1_MUT', 'PIK3CA_MUT', 'KRAS_MUT', 'ERBB2_MUT'
false_neg_ras = false_neg_maf_df.query('Hugo_Symbol in ["KRAS", "AKT1", "ERBB2", "PIK3CA"]')
tru_pos_ras = true_pos_maf.query('Hugo_Symbol in ["KRAS", "AKT1", "ERBB2", "PIK3CA"]')
# Remove duplicate cell-lines. Assume 1 cosmic variant supercedes lack of cosmic
false_neg_dup = false_neg_ras.groupby('Tumor_Sample_Barcode')['isCOSMIChotspot']

cosmic_false_neg = (
    false_neg_dup.value_counts()
    .reset_index(name='count')
    .sort_values(by='isCOSMIChotspot')
    .drop_duplicates(subset='Tumor_Sample_Barcode', keep='last')
)

n_false_neg_cosmic = cosmic_false_neg['isCOSMIChotspot'].value_counts()[1]
n_false_neg_obs = cosmic_false_neg['isCOSMIChotspot'].shape[0]
# What is the proportion of COSMIC variant (True) to non-COSMIC variant (False)
print(cosmic_false_neg['isCOSMIChotspot'].value_counts())
cosmic_false_neg['isCOSMIChotspot'].value_counts(normalize=True)
# Remove duplicate cell-lines. Assume 1 cosmic variant supercedes lack of cosmic
tru_pos_dup = tru_pos_ras.groupby('Tumor_Sample_Barcode')['isCOSMIChotspot']

cosmic_true_pos = (
    tru_pos_dup.value_counts()
    .reset_index(name='count')
    .sort_values(by='isCOSMIChotspot')
    .drop_duplicates(subset='Tumor_Sample_Barcode', keep='last')
)

n_true_pos_cosmic = cosmic_true_pos['isCOSMIChotspot'].value_counts()[1]
n_tru_pos_obs = cosmic_true_pos['isCOSMIChotspot'].shape[0]
print(cosmic_true_pos['isCOSMIChotspot'].value_counts())
cosmic_true_pos['isCOSMIChotspot'].value_counts(normalize=True)
# Test if the proportions of observed COSMIC variants are significantly different
# between true positives and false negatives. All of these samples have Ras mutations.
# The question is asking if the proportion of Ras variants annotated in the COSMIC
# database is lower in False negative tumors (those the classifier predicted as
# Ras wild-type) than in True positive tumors
cosmic_prop_chi = proportions_chisquare(count = [n_false_neg_cosmic, n_true_pos_cosmic],
                                        nobs = [n_false_neg_obs, n_tru_pos_obs])

print('Chi Square = {}, p value = {}'.format(cosmic_prop_chi[0], cosmic_prop_chi[1]))
print('There is a no significant difference in the expected proportion of COSMIC variants.')<jupyter_output>Chi Square = 1.7817635604283832, p value = 0.1819329204145474
There is a no significant difference in the expected proportion of COSMIC variants.<jupyter_text>Add CCLE Variant Scores (nucleotide and amino acid) to Supplementary Data Files<jupyter_code># Load TCGA PanCanAtlas Core Ras Pathway genes
pi3k_genes_file = os.path.join('..', 'data', 'pi3k_genes.csv')
pi3k_core_df = pd.read_csv(pi3k_genes_file)
pi3k_core_df.head()
# Subset MAF file to Ras pathway variants and merge with CCLE classifier scores
pi3k_pathway_genes = pi3k_core_df['genes'].tolist()
all_common_lines = set(ccle_maf_df.index).intersection(set(ccle_full_df.index))

# Subset to common cell lines
subset_maf = ccle_maf_df.loc[list(all_common_lines), :]
subset_maf = (
    subset_maf.query('Hugo_Symbol in @pi3k_pathway_genes')
    .loc[:, ['Hugo_Symbol', 'Protein_Change', 'cDNA_Change']]
    .merge(ccle_full_df, left_index=True, right_index=True)
)

subset_maf.head(3)
# Get the mean classifier scores for CCLE nucleotide variants
mean_nuc_data = (
    pd.DataFrame(subset_maf
                 .groupby(['cDNA_Change', 'Hugo_Symbol'])['weight']
                 .mean())
)
mean_nuc_data.columns = ['ccle_mean_weight']
mean_nuc_data = mean_nuc_data.reset_index()

# Get the sd classifier scores for CCLE variants
sd_nuc_data = (
    pd.DataFrame(subset_maf
                 .groupby(['cDNA_Change', 'Hugo_Symbol'])['weight']
                 .std())
)
sd_nuc_data.columns = ['ccle_sd_weight']
sd_nuc_data = sd_nuc_data.reset_index()

# Counts of CCLE variants altering amino acids
count_nuc_data = (
    pd.DataFrame(subset_maf
                 .groupby(['cDNA_Change', 'Hugo_Symbol'])['weight']
                 .count())
)
count_nuc_data.columns = ['ccle_count']
count_nuc_data = count_nuc_data.reset_index()
# Merge protein data
nuc_merge_on = ['Hugo_Symbol', 'cDNA_Change']
nuc_change_df = (
    mean_nuc_data.merge(sd_nuc_data,
                        left_on=nuc_merge_on, right_on=nuc_merge_on)
    .merge(count_nuc_data, left_on=nuc_merge_on, right_on=nuc_merge_on)
)

nuc_change_df.sort_values('ccle_count').tail(5)
data_s4_file = os.path.join('..', 'classifiers', 'ERBB2_PIK3CA_KRAS_AKT1', 'tables',
                            'nucleotide_mutation_scores.tsv')
data_s4_df = pd.read_table(data_s4_file)
# Merge the CCLE nucleotide scores
data_s4_df = data_s4_df.merge(nuc_change_df, left_on = ['Hugo_Symbol', 'HGVSc'],
                                    right_on = ['Hugo_Symbol', 'cDNA_Change'],
                              how='outer')
updated_data_s4_df = data_s4_df.sort_values(by='count', ascending=False)
updated_data_s4_df.head()
updated_data_s4_file = os.path.join('..', 'classifiers', 'ERBB2_PIK3CA_KRAS_AKT1', 'tables', 'updated_Data_S4.csv')
updated_data_s4_df.to_csv(updated_data_s4_file, sep=',', index=False)
# Get the mean classifier scores for CCLE variants
mean_protein_data = (
    pd.DataFrame(subset_maf
                 .groupby(['Protein_Change', 'Hugo_Symbol'])['weight']
                 .mean())
)
mean_protein_data.columns = ['ccle_mean_weight']
mean_protein_data = mean_protein_data.reset_index()

# Get the sd classifier scores for CCLE variants
sd_protein_data = (
    pd.DataFrame(subset_maf
                 .groupby(['Protein_Change', 'Hugo_Symbol'])['weight']
                 .std())
)
sd_protein_data.columns = ['ccle_sd_weight']
sd_protein_data = sd_protein_data.reset_index()

# Counts of CCLE variants altering amino acids
count_protein_data = (
    pd.DataFrame(subset_maf
                 .groupby(['Protein_Change', 'Hugo_Symbol'])['weight']
                 .count())
)
count_protein_data.columns = ['ccle_count']
count_protein_data = count_protein_data.reset_index()
# Merge protein data
merge_on = ['Hugo_Symbol', 'Protein_Change']
protein_change_df = (
    mean_protein_data.merge(sd_protein_data,
                            left_on=merge_on, right_on=merge_on)
    .merge(count_protein_data, left_on=merge_on, right_on=merge_on)
)

protein_change_df.sort_values('ccle_count').tail(5)
# Convert amino acid to 3 letters
protein_convert = [''.join([aa[x] if x in aa.keys() else x for x in y]) 
                   for y in protein_change_df['Protein_Change']]

protein_change_df = protein_change_df.assign(conversion = protein_convert)
data_s5_file = os.path.join('..', 'classifiers', 'ERBB2_PIK3CA_KRAS_AKT1', 'tables',
                            'amino_acid_mutation_scores.tsv')
data_s5_df = pd.read_table(data_s5_file)

# Merge the CCLE protein scores
data_s5_df = data_s5_df.merge(protein_change_df, left_on = ['Hugo_Symbol', 'HGVSp'],
                                    right_on = ['Hugo_Symbol', 'conversion'],
                              how='outer')

# Sort by the total number of mutations observed
updated_data_s5_df = (
    data_s5_df.drop(['Protein_Change'], axis=1).sort_values(by='count', ascending=False)
)

updated_data_s5_file = os.path.join('..', 'classifiers', 'ERBB2_PIK3CA_KRAS_AKT1','tables', 'updated_Data_S5.csv')
updated_data_s5_df.to_csv(updated_data_s5_file, sep=',', index=False)
updated_data_s5_df.head()<jupyter_output><empty_output><jupyter_text>CCLE - Pharmacologic EfficacyHere, we process drug efficacy data on the CCLE dataset. Data obtained from https://portals.broadinstitute.org/ccle/data (Pharmacologic profiling) (signin required).A processed `.tsv` file is output to be visualized in `scripts/viz/ras_ccle_pharmacology.R`.<jupyter_code># Load in pharmacological results
pharm_file = os.path.join('..', 'data', 'CCLE_NP24.2009_Drug_data_2015.02.24.csv')
pharm_df = pd.read_csv(pharm_file, index_col=0)
pharm_df = pharm_df.assign(tissue = [' '.join(x[1:]) for x in pharm_df.index.str.split('_')])

pharm_full_df = pharm_df.merge(ccle_full_df, left_index=True, right_index=True)
pharm_full_df.head()
ccle_full_df.head()
ccle_full_df.shape

common_celllines_pharm = set(ccle_full_df.index).intersection(set(pharm_df.index))
print('There are {} cell lines in common'.format(len(common_celllines_pharm)))

pharm_full_df['Compound'].value_counts()
pharm_full_df['tissue'].value_counts()
# What is the cell line tissue representation?
compound_heatmap = pd.pivot_table(pharm_full_df[['tissue', 'Compound']],
                                  columns='tissue', index='Compound',
                                  aggfunc=len)

compound_heatmap = pd.DataFrame(compound_heatmap.unstack()).reset_index()
compound_heatmap.columns = ['tissue', 'Compound', 'count']
compound_heatmap = compound_heatmap.sort_values(by=['tissue', 'Compound'])
gg.ggplot(compound_heatmap, gg.aes('factor(tissue)', 'factor(Compound)', fill='count')) + \
gg.geom_tile(gg.aes(width=.95, height=.95)) + \
gg.theme(axis_text_x=gg.element_text(rotation=90),
         panel_background=gg.element_rect(fill='white'))
# Write out pharm_full_df to file to plot in ggplot2
# plotnine does not include all the functionality required to create the plot
pharm_file = os.path.join('..', 'data', 'pharmacology_predictions_ERBB2_PIK3CA_KRAS_AKT1_ccle.tsv')
pharm_full_df.to_csv(pharm_file, sep='\t')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12395061728395061, 'fraction_numerical': 0.021459034792368124, 'mean_word_length': 4.825313807531381, 'pattern_counts': {'":': 3, '<': 31, '<?xml version=': 0, '>': 33, 'https://': 5, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 17, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '36858', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 8693, 'n_tokens_neox': 8102, 'n_words': 1570, 'path': 'scripts/pi3k_cell_line_predictions.ipynb', 'repo_name': 'nvk747/pancancer', 'size': '22275'}