From 7aafa72e10baf0427e543a664658c2d5cbae3333 Mon Sep 17 00:00:00 2001 From: Tanushree Tunstall Date: Tue, 14 Apr 2020 11:30:36 +0100 Subject: [PATCH] defined method for formatting mcsm_results --- mcsm/format_results.py | 421 ++++++++++++++++++++++------------------- mcsm/mcsm_results.py | 30 ++- 2 files changed, 245 insertions(+), 206 deletions(-) diff --git a/mcsm/format_results.py b/mcsm/format_results.py index d0c8e52..6670022 100755 --- a/mcsm/format_results.py +++ b/mcsm/format_results.py @@ -6,9 +6,9 @@ import os,sys import subprocess import argparse -import requests +#import requests import re -import time +#import time import pandas as pd from pandas.api.types import is_string_dtype from pandas.api.types import is_numeric_dtype @@ -43,7 +43,7 @@ datadir = homedir + '/' + 'git/Data' #======= # 1) result_urls (from outdir) outdir = datadir + '/' + drug + '/' + 'output' -in_filename = gene.lower() + '_mcsm_output.csv' #(outfile, from mcsm_results) +in_filename = gene.lower() + '_mcsm_output.csv' #(outfile, from mcsm_results.py) infile = outdir + '/' + in_filename print('Input filename:', in_filename , '\nInput path(from output dir):', outdir @@ -53,215 +53,256 @@ print('Input filename:', in_filename # output #======= outdir = datadir + '/' + drug + '/' + 'output' -out_filename = gene.lower() + '_complex_mcsm_norm.csv' +out_filename = gene.lower() + '_complex_mcsm_results.csv' outfile = outdir + '/' + out_filename print('Output filename:', out_filename , '\nOutput path:', outdir , '\n=============================================================') #======================================================================= -print('Reading input file') -mcsm_data = pd.read_csv(infile, sep = ',') +def format_mcsm_output(mcsm_outputcsv): + """ + @param mcsm_outputcsv: file containing mcsm results for all muts + (outfile from from mcsm_results.py) + @type string + + @return formatted mcsm output + @type pandas df + + """ + ############# + # Read file + ############# + mcsm_data = pd.read_csv(infile, sep = ',') + dforig_shape = mcsm_data.shape + print('dim of infile:', dforig_shape) + + ############# + # rename cols + ############# + # format colnames: all lowercase, remove spaces and use '_' to join + print('Assigning meaningful colnames i.e without spaces and hyphen and reflecting units' + , '\n===================================================================') + my_colnames_dict = {'Predicted Affinity Change': 'PredAffLog' # relevant info from this col will be extracted and the column discarded + , 'Mutation information': 'mutation_information' # {wild_type}{mutant_type} + , 'Wild-type': 'wild_type' # one letter amino acid code + , 'Position': 'position' # number + , 'Mutant-type': 'mutant_type' # one letter amino acid code + , 'Chain': 'chain' # single letter (caps) + , 'Ligand ID': 'ligand_id' # 3-letter code + , 'Distance to ligand': 'ligand_distance' # angstroms + , 'DUET stability change': 'duet_stability_change'} # in kcal/mol -mcsm_data.columns -# PredAffLog = affinity_change_log -# "DUETStability_Kcalpermol = DUET_change_kcalpermol -dforig_shape = mcsm_data.shape -print('dim of infile:', dforig_shape) + mcsm_data.rename(columns = my_colnames_dict, inplace = True) + #%%=========================================================================== + ################################# + # populate mutation_information + # col which is currently blank + ################################# + # populate mutation_information column:mcsm style muts {WT}{MUT} + print('Populating column : mutation_information which is currently empty\n', mcsm_data['mutation_information']) + mcsm_data['mutation_information'] = mcsm_data['wild_type'] + mcsm_data['position'].astype(str) + mcsm_data['mutant_type'] + print('checking after populating:\n', mcsm_data['mutation_information'] + , '\n===================================================================') -# change colnames to reflect units and no spaces, and replace '-' with '-' -print('Assigning meaningful colnames i.e without spaces and hyphen and reflecting units' - , '\n===================================================================') -my_colnames_dict = {'Predicted Affinity Change': 'PredAffLog' - , 'Mutation information': 'Mutationinformation' - , 'Wild-type': 'Wild_type' - , 'Position': 'Position' - , 'Mutant-type': 'Mutant_type' - , 'Chain': 'Chain' - , 'Ligand ID': 'LigandID' - , 'Distance to ligand': 'Dis_lig_Ang' - , 'DUET stability change': 'DUET_change_kcalpermol'} + # Remove spaces b/w pasted columns + print('removing white space within column: \mutation_information') + mcsm_data['mutation_information'] = mcsm_data['mutation_information'].str.replace(' ', '') + print('Correctly formatted column: mutation_information\n', mcsm_data['mutation_information'] + , '\n===================================================================') + #%%=========================================================================== + ############# + # sanity check: drop dupliate muts + ############# + # shouldn't exist as this should be eliminated at the time of running mcsm + print('Sanity check:' + , '\nChecking duplicate mutations') + if mcsm_data['mutation_information'].duplicated().sum() == 0: + print('PASS: No duplicate mutations detected (as expected)' + , '\nDim of data:', mcsm_data.shape + , '\n===============================================================') + else: + print('FAIL (but not fatal): Duplicate mutations detected' + , '\nDim of df with duplicates:', mcsm_data.shape + , 'Removing duplicate entries') + mcsm_data = mcsm_data.drop_duplicates(['mutation_information']) + print('Dim of data after removing duplicate muts:', mcsm_data.shape + , '\n===============================================================') + #%%=========================================================================== + ############# + # Create col: duet_outcome + ############# + # classification based on DUET stability values + print('Assigning col: duet_outcome based on DUET stability values') + print('Sanity check:') + # count positive values in the DUET column + c = mcsm_data[mcsm_data['duet_stability_change']>=0].count() + DUET_pos = c.get(key = 'duet_stability_change') + # Assign category based on sign (+ve : Stabilising, -ve: Destabilising, Mind the spelling (British spelling)) + mcsm_data['duet_outcome'] = np.where(mcsm_data['duet_stability_change']>=0, 'Stabilising', 'Destabilising') + mcsm_data['duet_outcome'].value_counts() + if DUET_pos == mcsm_data['duet_outcome'].value_counts()['Stabilising']: + print('PASS: DUET outcome assigned correctly') + else: + print('FAIL: DUET outcome assigned incorrectly' + , '\nExpected no. of stabilising mutations:', DUET_pos + , '\nGot no. of stabilising mutations', mcsm_data['duet_outcome'].value_counts()['Stabilising'] + , '\n===============================================================') + #%%=========================================================================== + ############# + # Extract numeric + # part of ligand_distance col + ############# + # Extract only the numeric part from col: ligand_distance + # number: '-?\d+\.?\d*' + mcsm_data['ligand_distance'] + print('extracting numeric part of col: ligand_distance') + mcsm_data['ligand_distance'] = mcsm_data['ligand_distance'].str.extract('(\d+\.?\d*)') + mcsm_data['ligand_distance'] + #%%=========================================================================== + ############# + # Create 2 columns: + # ligand_affinity_change and ligand_outcome + ############# + # the numerical and categorical parts need to be extracted from column: PredAffLog + # regex used + # numerical part: '-?\d+\.?\d*' + # categorocal part: '\b(\w+ing)\b' + print('Extracting numerical and categorical parts from the col: PredAffLog') + print('to create two columns: ligand_affinity_change and ligand_outcome' + , '\n===================================================================') -mcsm_data.rename(columns = my_colnames_dict, inplace = True) -mcsm_data.columns -#%%=========================================================================== -# populate mutationinformation column:mcsm style muts {WT}{MUT} -print('Populating column : Mutationinformation which is currently empty\n', mcsm_data['Mutationinformation']) -mcsm_data['Mutationinformation'] = mcsm_data['Wild_type'] + mcsm_data['Position'].astype(str) + mcsm_data['Mutant_type'] -print('checking after populating:\n', mcsm_data['Mutationinformation'] - , '\n===================================================================') + # 1) Extracting the predicted affinity change (numerical part) + mcsm_data['ligand_affinity_change'] = mcsm_data['PredAffLog'].str.extract('(-?\d+\.?\d*)', expand = True) + print(mcsm_data['ligand_affinity_change']) + + # 2) Extracting the categorical part (Destabillizing and Stabilizing) using word boundary ('ing') + #aff_regex = re.compile(r'\b(\w+ing)\b') + mcsm_data['ligand_outcome']= mcsm_data['PredAffLog'].str.extract(r'(\b\w+ing\b)', expand = True) + print(mcsm_data['ligand_outcome']) + print(mcsm_data['ligand_outcome'].value_counts()) -# Remove spaces b/w pasted columns -print('removing white space within column: \Mutationinformation') -mcsm_data['Mutationinformation'] = mcsm_data['Mutationinformation'].str.replace(' ', '') -print('Correctly formatted column: Mutationinformation\n', mcsm_data['Mutationinformation'] - , '\n===================================================================') -#%%=========================================================================== -# very important -print('Sanity check:' - , '\nChecking duplicate mutations') -if mcsm_data['Mutationinformation'].duplicated().sum() == 0: - print('PASS: No duplicate mutations detected (as expected)' - , '\nDim of data:', mcsm_data.shape - , '\n===============================================================') -else: - print('FAIL (but not fatal): Duplicate mutations detected' - , '\nDim of df with duplicates:', mcsm_data.shape - , 'Removing duplicate entries') - mcsm_data = mcsm_data.drop_duplicates(['Mutationinformation']) - print('Dim of data after removing duplicate muts:', mcsm_data.shape - , '\n===============================================================') -#%%=========================================================================== -# create DUET_outcome column: classification based on DUET stability values -print('Assigning col: DUET_outcome based on DUET stability values') -print('Sanity check:') -# count positive values in the DUET column -c = mcsm_data[mcsm_data['DUET_change_kcalpermol']>=0].count() -DUET_pos = c.get(key = 'DUET_change_kcalpermol') -# Assign category based on sign (+ve : Stabilising, -ve: Destabilising, Mind the spelling (British spelling)) -mcsm_data['DUET_outcome'] = np.where(mcsm_data['DUET_change_kcalpermol']>=0, 'Stabilising', 'Destabilising') -mcsm_data['DUET_outcome'].value_counts() -if DUET_pos == mcsm_data['DUET_outcome'].value_counts()['Stabilising']: - print('PASS: DUET outcome assigned correctly') -else: - print('FAIL: DUET outcome assigned incorrectly' - , '\nExpected no. of stabilising mutations:', DUET_pos - , '\nGot no. of stabilising mutations', mcsm_data['DUET_outcome'].value_counts()['Stabilising'] - , '\n===============================================================') -#%%=========================================================================== -# Extract only the numeric part from col: Dis_lig_Ang -# number: '-?\d+\.?\d*' -mcsm_data['Dis_lig_Ang'] -print('extracting numeric part of col: Dis_lig_Ang') -mcsm_data['Dis_lig_Ang'] = mcsm_data['Dis_lig_Ang'].str.extract('(\d+\.?\d*)') -mcsm_data['Dis_lig_Ang'] + ############# + # changing spelling: British + ############# + # ensuring spellings are consistent + american_spl = mcsm_data['ligand_outcome'].value_counts() + print('Changing to Bristish spellings for col: ligand_outcome') + mcsm_data['ligand_outcome'].replace({'Destabilizing': 'Destabilising', 'Stabilizing': 'Stabilising'}, inplace = True) + print(mcsm_data['ligand_outcome'].value_counts()) + british_spl = mcsm_data['ligand_outcome'].value_counts() + # compare series values since index will differ from spelling change + check = american_spl.values == british_spl.values + if check.all(): + print('PASS: spelling change successfull' + , '\nNo. of predicted affinity changes:\n', british_spl + , '\n===============================================================') + else: + print('FAIL: spelling change unsucessfull' + , '\nExpected:\n', american_spl + , '\nGot:\n', british_spl + , '\n===============================================================') + #%%=========================================================================== + ############# + # ensuring corrrect dtype columns + ############# + # check dtype in cols + print('Checking dtypes in all columns:\n', mcsm_data.dtypes + , '\n===================================================================') + print('Converting the following cols to numeric:' + , '\nligand_distance' + , '\nduet_stability_change' + , '\nligand_affinity_change' + , '\n===================================================================') + + # using apply method to change stabilty and affinity values to numeric + numeric_cols = ['duet_stability_change', 'ligand_affinity_change', 'ligand_distance'] + mcsm_data[numeric_cols] = mcsm_data[numeric_cols].apply(pd.to_numeric) + # check dtype in cols + print('checking dtype after conversion') + cols_check = mcsm_data.select_dtypes(include='float64').columns.isin(numeric_cols) + if cols_check.all(): + print('PASS: dtypes for selected cols:', numeric_cols + , '\nchanged to numeric' + , '\n===============================================================') + else: + print('FAIL:dtype change to numeric for selected cols unsuccessful' + , '\n===============================================================') + print(mcsm_data.dtypes) + #%%=========================================================================== -# changing dtype to numeric -#if is_numeric_dtype(mcsm_data['Dis_lig_Ang']): -# print('Data type is already numeric, doing nothing') -#else: -# print('Changing dtype in col: Dis_lig_Ang to numeric since Distance should be numeric') -## FIXME: either do it here, or in the end for all the required cols at once -#%%=========================================================================== -# create Lig_outcome column: classification based on affinity change values -# the numerical and categorical parts need to be extracted from column: PredAffLog -# regex used -# number: '-?\d+\.?\d*' -# category: '\b(\w+ing)\b' -print('Extracting numerical and categorical parts from the col: PredAffLog') -print('to create two columns: affinity_change_log and Lig_outcome' - , '\n===================================================================') + ############# + # scale duet values + ############# + # Rescale values in DUET_change col b/w -1 and 1 so negative numbers + # stay neg and pos numbers stay positive + duet_min = mcsm_data['duet_stability_change'].min() + duet_max = mcsm_data['duet_stability_change'].max() -# Extracting the predicted affinity change (numerical part) -mcsm_data['affinity_change_log'] = mcsm_data['PredAffLog'].str.extract('(-?\d+\.?\d*)', expand = True) -print(mcsm_data['affinity_change_log']) -# Extracting the categorical part (Destabillizing and Stabilizing) using word boundary ('ing') -#aff_regex = re.compile(r'\b(\w+ing)\b') -mcsm_data['Lig_outcome']= mcsm_data['PredAffLog'].str.extract(r'(\b\w+ing\b)', expand = True) -print(mcsm_data['Lig_outcome']) -print(mcsm_data['Lig_outcome'].value_counts()) -american_spl = mcsm_data['Lig_outcome'].value_counts() -print('Changing to Bristish spellings for col: Lig_outcome') -mcsm_data['Lig_outcome'].replace({'Destabilizing': 'Destabilising', 'Stabilizing': 'Stabilising'}, inplace = True) -print(mcsm_data['Lig_outcome'].value_counts()) -british_spl = mcsm_data['Lig_outcome'].value_counts() + duet_scale = lambda x : x/abs(duet_min) if x < 0 else (x/duet_max if x >= 0 else 'failed') -# since series object will have different names on account of our spelling change -# use .equals -check = american_spl.values == british_spl.values -if check.all(): - print('PASS: spelling change successfull' - , '\nNo. of predicted affinity changes:\n', british_spl - , '\n===============================================================') -else: - print('FAIL: spelling change unsucessfull' - , '\nExpected:\n', american_spl - , '\nGot:\n', british_spl - , '\n===============================================================') -#%%=========================================================================== -# check dtype in cols -print('Checking dtypes in all columns:\n', mcsm_data.dtypes - , '\n===================================================================') -print('Converting the following cols to numeric:' - , '\nDis_lig_Ang' - , '\nDUET_change_kcalpermol' - , '\naffinity_change_log' - , '\n===================================================================') -# using apply method to change stabilty and affinity values to numeric -numeric_cols = ['DUET_change_kcalpermol', 'affinity_change_log', 'Dis_lig_Ang'] -mcsm_data[numeric_cols] = mcsm_data[numeric_cols].apply(pd.to_numeric) -# check dtype in cols -print('checking dtype after conversion') -cols_check = mcsm_data.select_dtypes(include='float64').columns.isin(numeric_cols) -if cols_check.all(): - print('PASS: dtypes for selected cols:', numeric_cols - , '\nchanged to numeric' - , '\n===============================================================') -else: - print('FAIL:dtype change to numeric for selected cols unsuccessful' - , '\n===============================================================') -#mcsm_data['Dis_lig_Ang', 'affinity_change_log'].apply(is_numeric_dtype(mcsm_data['Dis_lig_Ang', 'affinity_change_log'])) -print(mcsm_data.dtypes) + mcsm_data['duet_scaled'] = mcsm_data['duet_stability_change'].apply(duet_scale) + print('Raw duet scores:\n', mcsm_data['duet_stability_change'] + , '\n---------------------------------------------------------------' + , '\nScaled duet scores:\n', mcsm_data['duet_scaled']) + + #%%=========================================================================== + ############# + # scale affinity values + ############# + # rescale values in affinity change col b/w -1 and 1 so negative numbers + # stay neg and pos numbers stay positive + aff_min = mcsm_data['ligand_affinity_change'].min() + aff_max = mcsm_data['ligand_affinity_change'].max() -#%%=========================================================================== + aff_scale = lambda x : x/abs(aff_min) if x < 0 else (x/aff_max if x >= 0 else 'failed') -#%%=========================================================================== -# Normalise the DUET and affinity change cols -#converter = lambda x : x*2 if x < 10 else (x*3 if x < 20 else x) -duet_min = mcsm_data['DUET_change_kcalpermol'].min() -duet_max = mcsm_data['DUET_change_kcalpermol'].max() - -converter = lambda x : x/abs(duet_min) if x < 0 else (x/duet_max if x >= 0 else 'failed') - -mcsm_data['DUET_change_kcalpermol'] -mcsm_data['ratioDUET'] = mcsm_data['DUET_change_kcalpermol'].apply(converter) -mcsm_data['ratioDUET'] -#%%=========================================================================== -# Normalise the affinity change cols -aff_min = mcsm_data['affinity_change_log'].min() -aff_max = mcsm_data['affinity_change_log'].max() - -converter = lambda x : x/abs(aff_min) if x < 0 else (x/aff_max if x >= 0 else 'failed') -#converter(mcsm_data['affinity_change_log']) - -mcsm_data['affinity_change_log'] -mcsm_data['ratioPredAff'] = mcsm_data['affinity_change_log'].apply(converter) -mcsm_data['ratioPredAff'] -#============================================================================= -# Removing PredAff log column as it is not needed? -print('Removing col: PredAffLog since relevant info has been extracted from it') -mcsm_dataf = mcsm_data.drop(columns = ['PredAffLog']) -#%%=========================================================================== -expected_cols_toadd = 4 -dforig_len = dforig_shape[1] -expected_cols = dforig_len + expected_cols_toadd -if len(mcsm_dataf.columns) == expected_cols: - print('PASS: formatting successful' - , '\nformatted df has expected no. of cols:', expected_cols - , '\n---------------------------------------------------------------' - , '\ncolnames:', mcsm_dataf.columns - , '\n----------------------------------------------------------------' - , '\ndtypes in cols:', mcsm_dataf.dtypes - , '\n----------------------------------------------------------------' - , '\norig data shape:', dforig_shape - , '\nformatted df shape:', mcsm_dataf.shape - , '\n===============================================================') -else: - print('FAIL: something went wrong in formatting df' - , '\nExpected no. of cols:', expected_cols - , '\nGot no. of cols:', len(mcsm_dataf.columns) - , '\nCheck formatting' - , '\n===============================================================') + mcsm_data['affinity_scaled'] = mcsm_data['ligand_affinity_change'].apply(aff_scale) + print('Raw affinity scores:\n', mcsm_data['ligand_affinity_change'] + , '\n---------------------------------------------------------------' + , '\nScaled affinity scores:\n', mcsm_data['affinity_scaled']) + #============================================================================= + # Removing PredAff log column as it is not needed? + print('Removing col: PredAffLog since relevant info has been extracted from it') + mcsm_dataf = mcsm_data.drop(columns = ['PredAffLog']) + #%%=========================================================================== + ############# + # sanity check before writing file + ############# + expected_cols_toadd = 4 + dforig_len = dforig_shape[1] + expected_cols = dforig_len + expected_cols_toadd + if len(mcsm_dataf.columns) == expected_cols: + print('PASS: formatting successful' + , '\nformatted df has expected no. of cols:', expected_cols + , '\n---------------------------------------------------------------' + , '\ncolnames:', mcsm_dataf.columns + , '\n----------------------------------------------------------------' + , '\ndtypes in cols:', mcsm_dataf.dtypes + , '\n----------------------------------------------------------------' + , '\norig data shape:', dforig_shape + , '\nformatted df shape:', mcsm_dataf.shape + , '\n===============================================================') + else: + print('FAIL: something went wrong in formatting df' + , '\nExpected no. of cols:', expected_cols + , '\nGot no. of cols:', len(mcsm_dataf.columns) + , '\nCheck formatting' + , '\n===============================================================') + return mcsm_dataf #%%============================================================================ +# call function +mcsm_df_formatted = format_mcsm_output(infile) + # writing file print('Writing formatted df to csv') -mcsm_dataf.to_csv(outfile, index = False) +mcsm_df_formatted.to_csv(outfile, index = False) print('Finished writing file:' , '\nFilename:', out_filename , '\nPath:', outdir - , '\nExpected no. of rows:', len(mcsm_dataf) - , '\nExpected no. of cols:', len(mcsm_dataf.columns) + , '\nExpected no. of rows:', len(mcsm_df_formatted) + , '\nExpected no. of cols:', len(mcsm_df_formatted) , '\n=============================================================') #%% #End of script diff --git a/mcsm/mcsm_results.py b/mcsm/mcsm_results.py index 03de9db..e72be7d 100755 --- a/mcsm/mcsm_results.py +++ b/mcsm/mcsm_results.py @@ -58,13 +58,10 @@ print('Output filename:', out_filename , '\nOutput path:', outdir , '\n=============================================================') -#%% global variables -#HOST = "http://biosig.unimelb.edu.au" -#PREDICTION_URL = f"{HOST}/mcsm_lig/prediction" #======================================================================= def fetch_results(urltextfile): """ - Extract results data from the results page + Extract results data using the prediction url @params result_page of request_results() @type response object @@ -90,36 +87,37 @@ def fetch_results(urltextfile): def build_result_dict(web_result_raw): """ - Format web results which is inconveniently preformatted! + Build dict of mcsm output for a single mutation + Format web results which is preformatted to enable building result dict # preformatted string object: Problematic! - - # roughly bring these to the same format as the other + # make format consistent - @params web_result_raw directly from html parser extraction + @params web_result_raw: directly from html parser extraction @type string @returns result dict @type {} """ - - # remove blank lines from output + +# remove blank lines from output mytext = os.linesep.join([s for s in web_result_raw.splitlines() if s]) - # Predicted affintiy change and DUET stability change cols - # are are split over multiple lines and Mutation information is empty! +# affinity change and DUET stability change cols are are split over +# multiple lines and Mutation information is empty! mytext = mytext.replace('ange:\n', 'ange: ') - #print(mytext) +# print(mytext) - # initiliase result_dict +# initiliase result_dict result_dict = {} for line in mytext.split('\n'): fields = line.split(':') - #print(fields) +# print(fields) if len(fields) > 1: # since Mutaton information is empty dict_entry = dict([(x, y) for x, y in zip(fields[::2], fields[1::2])]) result_dict.update(dict_entry) - return result_dict + return result_dict + #======================================================================= #%% call function #request_results(infile_url)