much development

mcsm_na/format_results_mcsm_na.py

@@ -51,7 +51,7 @@ def format_mcsm_na_output(mcsm_na_output_tsv):
     print('Assigning meaningful colnames'
             , '\n=======================================================')
     my_colnames_dict = {'PDB_FILE': 'pdb_file' # relevant info from this col will be extracted and the column discarded
-        , 'CHAIN': 'chain' # {wild_type}<position>{mutant_type}
+        , 'CHAIN': 'chain' 
         , 'WILD_RES': 'wild_type' # one letter amino acid code
         , 'RES_POS': 'position' # number
         , 'MUT_RES': 'mutant_type' # one letter amino acid code
@@ -65,8 +65,8 @@ def format_mcsm_na_output(mcsm_na_output_tsv):
     #############
     # create mutationinformation column
     #############    
-    mcsm_na_data['mutationinformation'] = mcsm_na_data['wild_type'] + mcsm_na_data.position.map(str) + mcsm_na_data['mutant_type']
+    #mcsm_na_data['mutationinformation'] = mcsm_na_data['wild_type'] + mcsm_na_data.position.map(str) + mcsm_na_data['mutant_type']
-
+    mcsm_na_data['mutationinformation'] = mcsm_na_data.loc[:,'wild_type'] + mcsm_na_data.loc[:,'position'].astype(int).apply(str) + mcsm_na_data.loc[:,'mutant_type']
 #%%===================================================================== 
     #############
     # Create col: mcsm_na_outcome
@@ -132,4 +132,3 @@ def format_mcsm_na_output(mcsm_na_output_tsv):
                                 , 'pdb_file']]
     return(mcsm_na_dataf)
 #%%##################################################################### 
-

scripts/combining_dfs.py

@@ -34,6 +34,11 @@ Created on Tue Aug  6 12:56:03 2019
 # Output: single csv of all 8 dfs combined
 # useful link
 # https://stackoverflow.com/questions/23668427/pandas-three-way-joining-multiple-dataframes-on-columns
+
+#%% FIXME: let the script proceed even if files don't exist!
+# i.e example below
+# '/home/tanu/git/Data/ethambutol/output/dynamut_results/embb_complex_dynamut_norm.csv'
+
 #=======================================================================
 #%% load packages
 import sys, os
@@ -48,7 +53,7 @@ homedir = os.path.expanduser('~')
 
 # set working dir
 os.getcwd()
-os.chdir(homedir + '/git/LSHTM_analysis/scripts')
+#os.chdir(homedir + '/git/LSHTM_analysis/scripts')
 os.getcwd()
 
 # FIXME: local imports
@@ -109,46 +114,80 @@ if not outdir:
 #=======
 # input
 #=======
-#in_filename_mcsm = gene.lower() + '_complex_mcsm_norm.csv' 
+gene_list_normal = ["pnca", "katg", "rpob", "alr"]
-in_filename_mcsm     = gene.lower() + '_complex_mcsm_norm_SAM.csv' # gidb
+
-in_filename_foldx    = gene.lower() + '_foldx.csv'
+if gene.lower() == "gid":
-in_filename_deepddg  = gene.lower() + '_ni_deepddg.csv' # change to decent filename and put it in the correct dir
+    print("\nReading mCSM file for gene:", gene)
-in_filename_dssp     = gene.lower() + '_dssp.csv'
+    in_filename_mcsm = gene.lower() + '_complex_mcsm_norm_SAM.csv'
-in_filename_kd       = gene.lower() + '_kd.csv'
+if gene.lower() == "embb":
-in_filename_rd       = gene.lower() + '_rd.csv'
+    print("\nReading mCSM file for gene:", gene)
-#in_filename_snpinfo = 'ns' + gene.lower() + '_snp_info_f.csv' # gwas f info
+    in_filename_mcsm = gene.lower() + '_complex_mcsm_norm1.csv'
-in_filename_afor     = gene.lower() + '_af_or.csv'
+if gene.lower() in gene_list_normal:
-#in_filename_afor_kin = gene.lower() + '_af_or_kinship.csv'
+   print("\nReading mCSM file for gene:", gene)
-infilename_dynamut   = gene.lower() + '_complex_dynamut_norm.csv'
+   in_filename_mcsm = gene.lower() + '_complex_mcsm_norm.csv' 
-infilename_dynamut2  = gene.lower() + '_complex_dynamut2_norm.csv'
-infilename_mcsm_na   = gene.lower() + '_complex_mcsm_na_norm.csv'
-infilename_mcsm_f_snps = gene.lower() + '_mcsm_formatted_snps.csv'
     
 infile_mcsm           = outdir + in_filename_mcsm
-infile_foldx   = outdir + in_filename_foldx
-infile_deepddg = outdir + in_filename_deepddg
-infile_dssp    = outdir + in_filename_dssp
-infile_kd      = outdir + in_filename_kd
-infile_rd      = outdir + in_filename_rd
-#infile_snpinfo = outdir + in_filename_snpinfo 
-infile_afor    = outdir + in_filename_afor
-#infile_afor_kin = outdir + in_filename_afor_kin
-infile_dynamut = outdir + 'dynamut_results/' + infilename_dynamut
-infile_dynamut2 = outdir + 'dynamut_results/dynamut2/' + infilename_dynamut2
-infile_mcsm_na = outdir + 'mcsm_na_results/' + infilename_mcsm_na
-infile_mcsm_f_snps = outdir + infilename_mcsm_f_snps
-
-# read csv
 mcsm_df               = pd.read_csv(infile_mcsm, sep = ',')
+
+in_filename_foldx     = gene.lower() + '_foldx.csv'
+infile_foldx          = outdir + in_filename_foldx
 foldx_df              = pd.read_csv(infile_foldx , sep = ',')
+
+in_filename_deepddg   = gene.lower() + '_ni_deepddg.csv' # change to decent filename and put it in the correct dir
+infile_deepddg        = outdir + in_filename_deepddg
 deepddg_df            = pd.read_csv(infile_deepddg, sep = ',')
+
+in_filename_dssp      = gene.lower() + '_dssp.csv'
+infile_dssp           = outdir + in_filename_dssp
 dssp_df               = pd.read_csv(infile_dssp, sep = ',')
+
+in_filename_kd        = gene.lower() + '_kd.csv'
+infile_kd             = outdir + in_filename_kd
 kd_df                 = pd.read_csv(infile_kd, sep = ',')
+
+in_filename_rd        = gene.lower() + '_rd.csv'
+infile_rd             = outdir + in_filename_rd
 rd_df                 = pd.read_csv(infile_rd, sep = ',')
+
+#in_filename_snpinfo   = 'ns' + gene.lower() + '_snp_info_f.csv' # gwas f info
+#infile_snpinfo        = outdir + in_filename_snpinfo 
+
+in_filename_afor      = gene.lower() + '_af_or.csv'
+infile_afor           = outdir + in_filename_afor
 afor_df               = pd.read_csv(infile_afor, sep = ',') 
-dynamut_df   = pd.read_csv(infile_dynamut, sep = ',')
+
+#in_filename_afor_kin  = gene.lower() + '_af_or_kinship.csv'
+#infile_afor_kin       = outdir + in_filename_afor_kin
+
+infilename_dynamut2    = gene.lower() + '_dynamut2_norm.csv'
+infile_dynamut2        = outdir + 'dynamut_results/dynamut2/' + infilename_dynamut2
 dynamut2_df            = pd.read_csv(infile_dynamut2, sep = ',')
-mcsm_na_df   = pd.read_csv(infile_mcsm_na, sep = ',')
+
+#------------------------------------------------------------
+# ONLY:for gene pnca and gid: End logic should pick this up!
+geneL_dy_na = ["pnca", "gid"]
+#if gene.lower() == "pnca" or "gid" :
+if gene.lower() in geneL_dy_na :
+    print("\nGene:", gene.lower()
+          , "\nReading Dynamut and mCSM_na files")
+    infilename_dynamut    = gene.lower() + '_dynamut_norm.csv' # gid
+    infile_dynamut        = outdir + 'dynamut_results/' + infilename_dynamut
+    dynamut_df            = pd.read_csv(infile_dynamut, sep = ',')
+    
+    infilename_mcsm_na    = gene.lower() + '_complex_mcsm_na_norm.csv' # gid
+    infile_mcsm_na        = outdir + 'mcsm_na_results/' + infilename_mcsm_na 
+    mcsm_na_df            = pd.read_csv(infile_mcsm_na, sep = ',')
+
+# ONLY:for gene embb and alr: End logic should pick this up!
+geneL_ppi2 = ["embb", "alr"]
+#if gene.lower() == "embb" or "alr":
+if gene.lower() in "embb" or "alr":
+    infilename_mcsm_ppi2   = gene.lower() + '_complex_mcsm_ppi2_norm.csv'
+    infile_mcsm_ppi2       = outdir + 'mcsm_ppi2/' + infilename_mcsm_ppi2
+    mcsm_ppi2_df           = pd.read_csv(infile_mcsm_ppi2, sep = ',')
+#--------------------------------------------------------------
+infilename_mcsm_f_snps = gene.lower() + '_mcsm_formatted_snps.csv'
+infile_mcsm_f_snps     = outdir + infilename_mcsm_f_snps
 mcsm_f_snps            = pd.read_csv(infile_mcsm_f_snps, sep = ',', names = ['mutationinformation'], header = None)
 
 #=======
@@ -158,12 +197,6 @@ out_filename_comb = gene.lower() + '_all_params.csv'
 outfile_comb =  outdir + out_filename_comb
 print('Output filename:', outfile_comb
       , '\n===================================================================')
-
-o_join = 'outer'
-l_join = 'left'
-r_join = 'right'
-i_join = 'inner'
-
 # end of variable assignment for input and output files
 #%%############################################################################  
 #=====================
@@ -292,6 +325,44 @@ else:
           , '\n======================================================')
     sys.exit()
 
+#--------------------------
+# check if >1 chain
+#--------------------------
+deepddg_df.loc[:,'chain_id'].value_counts()
+
+if len(deepddg_df.loc[:,'chain_id'].value_counts()) > 1:
+    print("\nChains detected: >1"
+          , "\nGene:", gene
+          , "\nChains:", deepddg_df.loc[:,'chain_id'].value_counts().index)
+    
+#--------------------------
+# subset chain
+#--------------------------
+if gene.lower() == "embb":
+    sel_chain = "B"
+else:
+    sel_chain = "A"
+    
+deepddg_df = deepddg_df[deepddg_df['chain_id'] == sel_chain]
+    
+#--------------------------
+# Check for duplicates
+#--------------------------
+if len(deepddg_df['mutationinformation'].duplicated().value_counts())> 1:
+    print("\nFAIL: Duplicates detected in DeepDDG infile"
+          , "\nNo. of duplicates:"
+          , deepddg_df['mutationinformation'].duplicated().value_counts()[1]
+          , "\nformat deepDDG infile before proceeding")
+    sys.exit()
+else:
+    print("\nPASS: No duplicates detected in DeepDDG infile")
+
+#--------------------------
+# Drop chain id col as other targets don't have itCheck for duplicates
+#--------------------------
+col_to_drop = ['chain_id']
+deepddg_df = deepddg_df.drop(col_to_drop, axis = 1)
+
 #%%=============================================================================
 # Now merges begin
 #%%=============================================================================
@@ -311,28 +382,83 @@ get_aa_3lower(df = mcsm_df
 #mcsm_df.columns = mcsm_df.columns.str.lower()
 # foldx_df.shape
 
-#mcsm_foldx_dfs = combine_dfs_with_checks(mcsm_df, foldx_df, my_join = o_join)
+#mcsm_foldx_dfs = combine_dfs_with_checks(mcsm_df, foldx_df, my_join = "outer")
 merging_cols_m1  = detect_common_cols(mcsm_df, foldx_df)
-mcsm_foldx_dfs = pd.merge(mcsm_df, foldx_df, on = merging_cols_m1,  how =  o_join)
+mcsm_foldx_dfs   = pd.merge(mcsm_df
+                          , foldx_df
+                          , on = merging_cols_m1
+                          , how =  "outer")
 ncols_m1 = len(mcsm_foldx_dfs.columns)
 
 print('\n\nResult of first merge:', mcsm_foldx_dfs.shape
       , '\n===================================================================')
 mcsm_foldx_dfs[merging_cols_m1].apply(len)
 mcsm_foldx_dfs[merging_cols_m1].apply(len) == len(mcsm_foldx_dfs)
+
+#%% for embB and any other targets where mCSM-lig hasn't run for 
+# get the empty cells to be full of meaningful info
+if mcsm_foldx_dfs.loc[:,'wild_type': 'mut_aa_3lower'].isnull().values.any():  
+    print ("NAs detected in mcsm cols after merge")
+    
+    ##############################
+    # Extract relevant col values 
+    # code to one
+    ##############################
+
+    # wt_reg = r'(^[A-Z]{1})'
+    # print('wild_type:', wt_reg)
+
+    # mut_reg = r'[0-9]+(\w{1})$'
+    # print('mut type:', mut_reg)
+    mcsm_foldx_dfs['wild_type']     = mcsm_foldx_dfs.loc[:,'mutationinformation'].str.extract(r'(^[A-Z]{1})')
+    mcsm_foldx_dfs['position']      = mcsm_foldx_dfs.loc[:,'mutationinformation'].str.extract(r'([0-9]+)')
+    mcsm_foldx_dfs['mutant_type']   = mcsm_foldx_dfs.loc[:,'mutationinformation'].str.extract(r'[0-9]+([A-Z]{1})$')
+    
+    # BEWARE: Bit of logic trap i.e if nan comes first
+    # in chain column, then nan will be populated!
+    #df['foo'] = df['chain'].unique()[0]
+    mcsm_foldx_dfs['chain'] = np.where(mcsm_foldx_dfs[['chain']].isnull().all(axis=1)
+                                        , mcsm_foldx_dfs['chain'].unique()[0]
+                                        , mcsm_foldx_dfs['chain'])
+    
+    mcsm_foldx_dfs['ligand_id'] = np.where(mcsm_foldx_dfs[['ligand_id']].isnull().all(axis=1)
+                                        , mcsm_foldx_dfs['ligand_id'].unique()[0]
+                                        , mcsm_foldx_dfs['ligand_id'])
+    #--------------------------------------------------------------------------
+    
+    mcsm_foldx_dfs['wild_pos']       = mcsm_foldx_dfs.loc[:,'wild_type'] + mcsm_foldx_dfs.loc[:,'position'].astype(int).apply(str)
+    mcsm_foldx_dfs['wild_chain_pos'] = mcsm_foldx_dfs.loc[:,'wild_type'] + mcsm_foldx_dfs.loc[:,'chain'] +  mcsm_foldx_dfs.loc[:,'position'].astype(int).apply(str)
+    
+    #############
+    # Map 1 letter 
+    # code to 3Upper
+    #############
+    # initialise a sub dict that is lookup dict for 
+    # 3-LETTER aa code to 1-LETTER aa code
+    lookup_dict = dict()
+    for k, v in oneletter_aa_dict.items():
+        lookup_dict[k] = v['three_letter_code_lower']
+        wt = mcsm_foldx_dfs['wild_type'].squeeze() # converts to a series that map works on
+        mcsm_foldx_dfs['wt_aa_3lower'] = wt.map(lookup_dict)   
+        mut = mcsm_foldx_dfs['mutant_type'].squeeze()
+        mcsm_foldx_dfs['mut_aa_3lower'] = mut.map(lookup_dict)
+    
 #%%
 print('==================================='
       , '\nSecond merge: mcsm_foldx_dfs + deepddg'
       , '\n===================================')
 
-#deepddg_df =  pd.read_csv(infile_deepddg, sep = ',')
-#deepddg_df.columns
-
 # merge with mcsm_foldx_dfs and deepddg_df
-mcsm_foldx_deepddg_dfs = pd.merge(mcsm_foldx_dfs, deepddg_df, on = 'mutationinformation',  how = l_join)
+mcsm_foldx_deepddg_dfs = pd.merge(mcsm_foldx_dfs
+                                  , deepddg_df
+                                  , on = 'mutationinformation'
+                                  , how = "left")
 mcsm_foldx_deepddg_dfs['deepddg_outcome'].value_counts()
 
 ncols_deepddg_merge = len(mcsm_foldx_deepddg_dfs.columns)
+
+mcsm_foldx_deepddg_dfs['position'] = mcsm_foldx_deepddg_dfs['position'].astype('int64')
+
 #%%============================================================================
 print('==================================='
       , '\Third merge: dssp + kd'
@@ -342,9 +468,12 @@ dssp_df.shape
 kd_df.shape
 rd_df.shape
 
-#dssp_kd_dfs = combine_dfs_with_checks(dssp_df, kd_df, my_join = o_join)
+#dssp_kd_dfs = combine_dfs_with_checks(dssp_df, kd_df, my_join = "outer")
 merging_cols_m2 = detect_common_cols(dssp_df, kd_df)
-dssp_kd_dfs = pd.merge(dssp_df, kd_df, on = merging_cols_m2,  how = o_join)
+dssp_kd_dfs = pd.merge(dssp_df
+                       , kd_df
+                       , on = merging_cols_m2
+                       , how = "outer")
 
 print('\n\nResult of third merge:', dssp_kd_dfs.shape
       , '\n===================================================================')
@@ -353,10 +482,12 @@ print('==================================='
       , '\nFourth merge: third merge + rd_df' 
       , '\ndssp_kd_dfs + rd_df'
       , '\n===================================')
-#dssp_kd_rd_dfs = combine_dfs_with_checks(dssp_kd_dfs, rd_df, my_join = o_join)
+#dssp_kd_rd_dfs = combine_dfs_with_checks(dssp_kd_dfs, rd_df, my_join = "outer")
 merging_cols_m3 = detect_common_cols(dssp_kd_dfs,  rd_df)
-dssp_kd_rd_dfs = pd.merge(dssp_kd_dfs, rd_df, on = merging_cols_m3
+dssp_kd_rd_dfs = pd.merge(dssp_kd_dfs
-                          , how = o_join)
+                          , rd_df
+                          , on = merging_cols_m3
+                          , how = "outer")
 
 ncols_m3 = len(dssp_kd_rd_dfs.columns)
 
@@ -369,23 +500,40 @@ print('======================================='
       , '\nFifth merge: Second merge + fourth merge'
       , '\nmcsm_foldx_dfs + dssp_kd_rd_dfs'
       , '\n=======================================')
-#combined_df = combine_dfs_with_checks(mcsm_foldx_dfs, dssp_kd_rd_dfs, my_join = i_join)
+
+#combined_df = combine_dfs_with_checks(mcsm_foldx_dfs, dssp_kd_rd_dfs, my_join = "inner")
 #merging_cols_m4 = detect_common_cols(mcsm_foldx_dfs, dssp_kd_rd_dfs)
-#combined_df = pd.merge(mcsm_foldx_dfs, dssp_kd_rd_dfs, on = merging_cols_m4, how  = i_join)
+#combined_df = pd.merge(mcsm_foldx_dfs, dssp_kd_rd_dfs, on = merging_cols_m4, how  = "inner")
 #combined_df_expected_cols = ncols_m1 + ncols_m3 - len(merging_cols_m4)
 
 # with deepddg values
 merging_cols_m4 = detect_common_cols(mcsm_foldx_deepddg_dfs, dssp_kd_rd_dfs)
-combined_df = pd.merge(mcsm_foldx_deepddg_dfs, dssp_kd_rd_dfs, on = merging_cols_m4, how  = i_join)
+combined_df = pd.merge(mcsm_foldx_deepddg_dfs
+                       , dssp_kd_rd_dfs
+                       , on = merging_cols_m4
+                       , how  = "inner")
 
 combined_df_expected_cols = ncols_deepddg_merge + ncols_m3 - len(merging_cols_m4)
 
-if len(combined_df) == len(mcsm_df) and len(combined_df.columns) == combined_df_expected_cols:
+# FIXME: check logic, doesn't effect anything else!
+if not gene == "embB":
+    print("\nGene is:", gene)
+    if len(combined_df) == len(mcsm_df) and len(combined_df.columns) == combined_df_expected_cols:
         print('PASS: successfully combined 5 dfs'
               , '\nNo. of rows combined_df:', len(combined_df)
               , '\nNo. of cols combined_df:', len(combined_df.columns))   
-else:
+    else:
-    sys.exit('FAIL: check individual df merges')
+        #sys.exit('FAIL: check individual df merges')
+        print("\nGene is:", gene
+              , "\ncombined_df length:", len(combined_df)
+              , "\nmcsm_df_length:", len(mcsm_df)
+              )
+        if len(combined_df.columns) == combined_df_expected_cols:
+            print('PASS: successfully combined 5 dfs'
+                  , '\nNo. of rows combined_df:', len(combined_df)
+                  , '\nNo. of cols combined_df:', len(combined_df.columns))
+        else:
+            sys.exit('FAIL: check individual merges')        
             
 print('\nResult of Fourth merge:', combined_df.shape
       , '\n===================================================================')
@@ -401,7 +549,7 @@ combined_df['chain'].equals(combined_df['chain_id'])
 combined_df['wild_type'].equals(combined_df['wild_type_kd']) # has nan
 combined_df['wild_type'].equals(combined_df['wild_type_dssp'])
 
-#sanity check
+# sanity check
 foo = combined_df[['wild_type', 'wild_type_kd', 'wt_3letter_caps', 'wt_aa_3lower', 'mut_aa_3lower']] 
 
 # Drop cols
@@ -455,7 +603,11 @@ afor_df = afor_df.drop(['position'], axis = 1)
 afor_cols = afor_df.columns
 
 # merge 
-combined_stab_afor = pd.merge(combined_df_clean, afor_df, on = merging_cols_m5, how  = l_join)
+combined_stab_afor = pd.merge(combined_df_clean
+                              , afor_df
+                              , on = merging_cols_m5
+                              , how  = "left")
+
 comb_afor_df_cols = combined_stab_afor.columns
 
 comb_afor_expected_cols = len(combined_df_clean.columns) + len(afor_df.columns) - len(merging_cols_m5)
@@ -467,15 +619,23 @@ if len(combined_stab_afor) == len(combined_df_clean) and len(combined_stab_afor.
 else:
     sys.exit('\nFAIL: check individual df merges')
 
-print('\n\nResult of Fourth merge:', combined_stab_afor.shape
+print('\n\nResult of Fifth merge:', combined_stab_afor.shape
       , '\n===================================================================')
 
 combined_stab_afor[merging_cols_m5].apply(len)
 combined_stab_afor[merging_cols_m5].apply(len) == len(combined_stab_afor)
 
-if len(combined_stab_afor) -  combined_stab_afor['mutation'].isna().sum() == len(afor_df):
+if (len(combined_stab_afor) - combined_stab_afor['mutation'].isna().sum()) == len(afor_df):
-     print('\nPASS: Merge successful for af and or'
+     print('\nPASS: Merge successful for af and or with matched numbers')
-          , '\nNo. of nsSNPs with valid ORs: ', len(afor_df))
+
+if len(combined_stab_afor) - combined_stab_afor['mutation'].isna().sum() == len(afor_df)-len(afor_df[~afor_df['mutation'].isin(combined_stab_afor['mutation'])]):
+    print("\nMismatched numbers, OR df has extra snps not found in mcsm df"
+          , "\nNo. of nsSNPs with valid ORs:", len(afor_df)
+          , "\nNo. of mcsm nsSNPs: ", len(combined_df_clean) 
+          , "\nNo. of OR nsSNPs not in mCSM df:"
+          , len(afor_df[~afor_df['mutation'].isin(combined_stab_afor['mutation'])])
+          , "\nWriting these mutations to file:")
+    orsnps_notmcsm = afor_df[~afor_df['mutation'].isin(combined_stab_afor['mutation'])]
 else:
     sys.exit('\nFAIL: merge unsuccessful for af and or')    
     
@@ -486,7 +646,7 @@ outfile_comb_afor =  outdir + '/' + out_filename_comb_afor
 print('Output filename:', outfile_comb_afor
       , '\n===================================================================')
 
-# # write csv
+# write csv
 print('Writing file: combined stability and afor')
 combined_stab_afor.to_csv(outfile_comb_afor, index = False)
 print('\nFinished writing file:'
@@ -494,7 +654,20 @@ print('\nFinished writing file:'
       , '\nNo. of cols:', combined_stab_afor.shape[1])
 #%%============================================================================
 # combine dynamut, dynamut2, and mcsm_na
-dfs_list = [dynamut_df, dynamut2_df, mcsm_na_df]
+#dfs_list = [dynamut_df, dynamut2_df, mcsm_na_df] # gid
+
+if gene.lower() == "pnca":
+    dfs_list = [dynamut_df, dynamut2_df]
+if gene.lower() == "gid":
+    dfs_list = [dynamut_df, dynamut2_df, mcsm_na_df]
+if gene.lower() == "embb":
+    dfs_list = [dynamut2_df, mcsm_ppi2_df]
+if gene.lower() == "katg":
+    dfs_list = [dynamut2_df]
+if gene.lower() == "rpob":
+    dfs_list = [dynamut2_df]
+if gene.lower() == "alr":
+    dfs_list = [dynamut2_df, mcsm_ppi2_df]
 
 dfs_merged = reduce(lambda  left,right: pd.merge(left
                                                 , right
@@ -514,7 +687,7 @@ len(combined_stab_afor.columns)
 combined_all_params = pd.merge(combined_stab_afor
                                , dfs_merged_clean
                                , on = merging_cols_m6
-                               , how  = i_join)
+                               , how  = "inner")
 
 expected_ncols = len(dfs_merged_clean.columns) + len(combined_stab_afor.columns) - len(merging_cols_m6)
 expected_nrows = len(combined_stab_afor)

scripts/data_extraction.py

@@ -70,7 +70,6 @@ arg_parser.add_argument('-m', '--make_dirs', help = 'Make dir for input and outp
 
 arg_parser.add_argument('--debug', action ='store_true', help = 'Debug Mode')
 
-
 args = arg_parser.parse_args()
 #=======================================================================
 #%% variable assignment: input and output paths & filenames

scripts/deepddg_format.py

@@ -117,12 +117,20 @@ deepddg_df['deepddg_outcome'].value_counts()
 len(deepddg_df.loc[deepddg_df['deepddg'] < 0])
 len(deepddg_df.loc[deepddg_df['deepddg'] >= 0])
                
+#----------------------------------------------
 # drop extra columns to allow clean merging
-deepddg_short_df = deepddg_df.drop(['chain_id', 'wild_type_deepddg', 'position', 'mutant_type_deepddg'], axis = 1)
+#----------------------------------------------
+#deepddg_short_df = deepddg_df.drop(['chain_id', 'wild_type_deepddg', 'position', 'mutant_type_deepddg'], axis = 1)
+
+#----------------------------------------------
+# embb (where gene-target has > 1 chain)
+# include chain else the numbering will be messed up!
+#----------------------------------------------
+deepddg_short_df = deepddg_df.drop(['wild_type_deepddg', 'position', 'mutant_type_deepddg'], axis = 1)
 
 # rearrange columns
 deepddg_short_df.columns
-deepddg_short_df = deepddg_short_df[["mutationinformation", "deepddg", "deepddg_outcome"]]
+deepddg_short_df = deepddg_short_df[["chain_id", "mutationinformation", "deepddg", "deepddg_outcome"]]
 
 #%% combine with mcsm snps 
 deepddg_mcsm_muts_dfs = pd.merge(deepddg_short_df

scripts/rd_df.py

@@ -45,8 +45,6 @@ arg_parser.add_argument('--debug', action='store_true', help = 'Debug Mode')
 args = arg_parser.parse_args()
 #=======================================================================
 #%% variable assignment: input and output 
-#drug = 'pyrazinamide'
-#gene = 'pncA'
 drug = args.drug
 gene = args.gene
 gene_match = gene + '_p.'