<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Permuation based Mutual Exclusive and Co-occurrence for group of items
Step1:  This is the python implementation of the MutEx algorithm described in the PhD Thesis "Analysis and visualization of multidimensional cancer genomics data" 
Step2:  We import mutex and create a MutEx instance
The MutEx has to be created with the full data matrix - so it has the correct background event (mutation) rate for both samples and genes
Step3:  Example run
Step4:  Example run , multi-threaded
Step5:  Example with many groups
Step6:  Big matrix
Step7:  Comparing simulation with 'original' data
We compute the mutational burden per sample and do 500 simulations of the mutation matrix
Step8:  The 500 simulations of the mutation matrix - parallel execution
Step9:  calculate mean and standard deviation of simulated mutational burden for each sample
Step10:  We plot the simulated mean with the standard deviation as well as the 'original' mutational burden 
We can see that the mean and the original mutational burden are overlapping very well with 500 simulated sets
<ASSISTANT_TASK:>
Python Code:
import numpy as np
import pandas as pd
import scipy.sparse as sparse
import random
Explanation: Permuation based Mutual Exclusive and Co-occurrence for group of items
End of explanation
row, col = 100, 100
np.random.seed(77)
df = pd.DataFrame(sparse.random(row, col, density=0.15).A).apply(np.ceil)
df.loc[0] = [1 if x < 20 else 0 for x in range(0, df.shape[1])]
df.loc[1] = [1 if x > 13 and x < 35 else 0 for x in range(0, df.shape[1])]
df.loc[2] = [1 if x > 80 else 0 for x in range(0, df.shape[1])]
df.columns = ['s' + str(x) for x in df.columns]
df.index = ['gene' + str(x) for x in df.index]
pd.set_option('display.max_columns', 1000)
df.head()
Explanation: This is the python implementation of the MutEx algorithm described in the PhD Thesis "Analysis and visualization of multidimensional cancer genomics data" : http://www.tdx.cat/handle/10803/301436
The same algorithm is also implement in Java in the Gitools software: http://www.gitools.org
Create a random sparse DataFrame reperesenting mutations
End of explanation
from mutex import MutEx
m = MutEx(background=df, permutations=2000)
Explanation: We import mutex and create a MutEx instance
The MutEx has to be created with the full data matrix - so it has the correct background event (mutation) rate for both samples and genes
End of explanation
%time r = m.calculate(['gene4', 'gene5', 'gene6'], parallel=False)
print(r)
Explanation: Example run
End of explanation
%time r = m.calculate(['gene0', 'gene1', 'gene2'])
print(r)
Explanation: Example run , multi-threaded
End of explanation
#some random groups and put it in a generator (alternatively list)
random.seed(18)
group_generator = (random.sample( df.index.tolist(), random.sample([2,3,4], 1)[0]) for x in range(10) )
result_list = [m.calculate(g) for g in group_generator]
result_df = pd.DataFrame.from_records([r.__dict__ for r in result_list])
result_df[['items','coverage', 'signal', 'mutex_pvalue', 'co_occurence_pvalue', 'signal_coverage_ratio', 'mean_sim_coverage']].sort_values('mutex_pvalue')
Explanation: Example with many groups
End of explanation
row, col = 10000, 800
np.random.seed(77)
df = pd.DataFrame(sparse.random(row, col).A).apply(np.ceil)
df.columns = ['s' + str(x) for x in df.columns]
df.index = ['gene' + str(x) for x in df.index]
pd.set_option('display.max_columns', 1000)
df.head()
m = MutEx(df, permutations=10000)
group = random.sample( df.index.tolist(), 20)
%time r = m.calculate(group)
r
Explanation: Big matrix
End of explanation
sample_mut_orig = df.apply(sum)
sample_mut_orig.head()
from pandas.tools.plotting import parallel_coordinates
import matplotlib.pylab as pylab
%matplotlib inline  
Explanation: Comparing simulation with 'original' data
We compute the mutational burden per sample and do 500 simulations of the mutation matrix
End of explanation
from functools import partial
import multiprocessing as mp
pool = mp.Pool(processes=mp.cpu_count())
def dosim(m,x):
    np.random.seed(x)
    return m._simulate_observations(df.apply(sum, axis=1)).apply(sum)
partial_simul = partial(dosim)
n = 500
%time simulated_results = pool.starmap(partial_simul.func, zip(itertools.repeat(m,n), range(n)) )
pool.close()  # we are not adding any more processes
pool.join()  # tell it to wait until all threads are done before going on
Explanation: The 500 simulations of the mutation matrix - parallel execution
End of explanation
sim_sample_muts = pd.concat(simulated_results, axis=1).transpose()
sim_sample_muts.columns = df.columns
sim_std = sim_sample_muts.apply(np.std)
sim_mean = sim_sample_muts.apply(np.mean)
sim_sample_muts.head()
Explanation: calculate mean and standard deviation of simulated mutational burden for each sample
End of explanation
sample_muts_stats = pd.concat([sim_mean, sim_mean + sim_std, sim_mean - sim_std, sample_mut_orig], axis=1).sort_values(0).transpose()
sample_muts_stats.head()
sample_muts_stats['simulated'] = ['mean', 'sim_mean + std', 'sim_mean - std', 'original']
pylab.rcParams['figure.figsize'] = 22, 10
parallel_coordinates(sample_muts_stats, 'simulated', axvlines=False)
Explanation: We plot the simulated mean with the standard deviation as well as the 'original' mutational burden 
We can see that the mean and the original mutational burden are overlapping very well with 500 simulated sets
End of explanation
<END_TASK>