<jupyter_start><jupyter_text>1. calculate $\Delta G$ of each microstate at pH 02. determine lowest energy tautomer of each charge state3. compare sequence of predicted dominant states to experimental states  1. Analyze microstates of Jaguar submission<jupyter_code>jaguartypei = SAMPL6DataProvider("typeI-raw-jaguar.csv", "typei", "Jaguar", bootstrap_options={"n_samples": 1})
exp = SAMPL6DataProvider("SAMPL6_experimental_pkas.csv", "exp", "Experiment", bootstrap_options={"n_samples": 3})
#"experimental_microstates.csv"
jaguartypei.__dict__
SM07_exp = exp.load("SM07")
SM07_jag = jaguartypei.load("SM07")

# show properties as example

SM07_jag.__dict__
microstate_data = []

for mol_ID in mol_IDs:
    jag_1mol = jaguartypei.load(mol_ID)
    microstate_IDs = jag_1mol.state_ids
    charges_of_microstates = jag_1mol.charges
    free_energies_of_microstates_pH0 = jag_1mol.free_energies[:,0]
    
    print(mol_ID)
    print("Microstates:\n", microstate_IDs)
    print("Charges:", charges_of_microstates)
    print("Free energies:\n",free_energies_of_microstates_pH0)
    print()
    
    for i, microstate_ID in enumerate(microstate_IDs):
    
        microstate_data.append({
            'Molecule ID': mol_ID,
            'Microstate ID': microstate_ID,
            'Charge': charges_of_microstates[i],
            '$\Delta$G (pH=0)': free_energies_of_microstates_pH0[i]
        })

df_microstate_data = pd.DataFrame(data=microstate_data)
df_microstate_data
df_microstate_data_charge0 = df_microstate_data[(df_microstate_data["Molecule ID"] == "SM07") & (df_microstate_data["Charge"] == 0)]
df_microstate_data_charge0
dominant_microstate = df_microstate_data_charge0.loc[df_microstate_data_charge0['$\Delta$G (pH=0)'].idxmin()]["Microstate ID"]
dominant_microstate
pred_mol_IDs = set(df_microstate_data["Molecule ID"].values)
pred_mol_IDs
# Create dataframe to store dominant predicted microstate ID of each charge

submission_ID = "jaguar"

df_dominant_ms = pd.DataFrame(columns=["receipt_id", "Molecule ID", "charge -4", "charge -3","charge -2", "charge -1","charge 0",
                                       "charge 1","charge 2", "charge 3", "charge 4"])
for i, mol_ID in enumerate(pred_mol_IDs): 
    df_dominant_ms.loc[i] = [submission_ID, mol_ID, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN ]
df_dominant_ms
# Populate dominant microstate dataframe

for mol_ID in pred_mol_IDs:
    df_1mol = df_microstate_data[df_microstate_data["Molecule ID"] == mol_ID]
    charges = set(df_1mol["Charge"].values)

    for charge in charges:
        df_1mol_1charge = df_1mol[df_1mol["Charge"] == charge]
        dominant_microstate = df_1mol_1charge.loc[df_1mol_1charge['$\Delta$G (pH=0)'].idxmin()]["Microstate ID"]
        df_dominant_ms.loc[(df_dominant_ms["Molecule ID"] == mol_ID), "charge {}".format(charge)] = dominant_microstate
        #print("Charge {}, {}".format(charge, dominant_microstate))

# Save as CSV
df_dominant_ms.to_csv("typeI_submission_dominant_microstate_collection_only_predictions.csv", index=False)
# Compare predicted dominant microstate collection to experimental dominant microstates
# Read predictied dominant microstate collection
df_pred_dom_ms = pd.read_csv("typeI_submission_dominant_microstate_collection_only_predictions.csv")
df_pred_dom_ms
# Read experimental dominant microstates
exp_microstates = pd.read_csv("experimental_microstates_with_charge.csv")
exp_microstates
# Organize experimental microstate in dominant microstate collection format

# Create empty dataframe
df_exp_dom_ms = pd.DataFrame(columns=[ "Molecule ID", "charge -4", "charge -3","charge -2", "charge -1","charge 0",
                                       "charge 1","charge 2", "charge 3", "charge 4"])
charges = np.arange(-4,4,1)
exp_mol_IDs = set(exp_microstates["Molecule ID"].values)

for i, mol_ID in enumerate(exp_mol_IDs): 
    df_exp_dom_ms.loc[i] = [mol_ID, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN ]

# Make molecule ID the index

df_exp_dom_ms = df_exp_dom_ms.set_index("Molecule ID", drop = False)
# Populate expermental dominant microstate collection
for i, mol_ID in enumerate(exp_mol_IDs): 
    for charge in charges:
        #print("mol {}, charge {}".format(mol_ID, charge))
        exp_ms_1mol = exp_microstates[exp_microstates["Molecule ID"] == mol_ID]
        #print("exp_ms_1mol:\n",exp_ms_1mol)
        
        # Check if charge exist in "Charge of A" column
        exp_charge_of_A = exp_ms_1mol["Charge of A"].values
        exp_charge_of_HA = exp_ms_1mol["Charge of HA"].values
        #print("exp_charge_of_A:", exp_charge_of_A)
        #print("exp_charge_of_HA:", exp_charge_of_HA)
        
        if charge in exp_charge_of_A:
            #print("Charges match A.")
            #exp_dom_ms_ID =  exp_ms_1mol["Microstate ID of A"].values
            exp_dom_ms_ID = exp_ms_1mol[exp_ms_1mol["Charge of A"] == charge]["Microstate ID of A"].values[0]
        elif charge in exp_charge_of_HA:
            #print("Charges match HA.")
            #exp_dom_ms_ID =  exp_ms_1mol["Microstate ID of HA"].values
            exp_dom_ms_ID = exp_ms_1mol[exp_ms_1mol["Charge of HA"] == charge]["Microstate ID of HA"].values[0]
        else:
            #print("Different charge")
            exp_dom_ms_ID = np.NaN
        
        #print("exp_dom_ms_ID:", exp_dom_ms_ID)
        
        df_exp_dom_ms.loc[mol_ID, "charge {}".format(charge)] = exp_dom_ms_ID

df_exp_dom_ms.to_csv("typeI_submission_dominant_microstate_collection_only_experiments.csv", index=False)
df_exp_dom_ms
# Create collection file with both experimental and predicted microstate IDs

# Read experimental data collection
df_exp_dom_ms = pd.read_csv("typeI_submission_dominant_microstate_collection_only_experiments.csv")
                           
# Read predicted microstates collections
df_pred_dom_ms = pd.read_csv("typeI_submission_dominant_microstate_collection_only_predictions.csv")

# TO-DO: Merge these dataframes create the collection file
df_pred_dom_ms
# Dominant microstate collection
df_dom_ms_col = pd.DataFrame(columns=["receipt_id","Molecule ID", "charge -4 (calc)", "charge -3 (calc)",
                                      "charge -2 (calc)", "charge -1 (calc)","charge 0 (calc)",
                                       "charge 1 (calc)","charge 2 (calc)", "charge 3 (calc)", "charge 4 (calc)",
                                     "charge -1 (exp)", "charge 0 (exp)", "charge 1 (exp)", "charge 2 (exp)"])
df_dom_ms_col
# Record data in a list to later convert to a dataframe
dominant_microstate_records = []

receipt_id_list = ["jaguar"]

for receipt_id in receipt_id_list:
    df_1method = df_pred_dom_ms[df_pred_dom_ms.receipt_id == receipt_id]
    for i, row in df_1method.iterrows():
        mol_ID = row["Molecule ID"]
        
        # Find the matching row from experimental dominant microstate dataframe using the mol_ID
        exp_row = df_exp_dom_ms[df_exp_dom_ms["Molecule ID"] == mol_ID]
        #print(exp_row)
        mol_ID_exp = exp_row["Molecule ID"].values[0]
        
        
        # Append to a list to later save as a CSV
        dominant_microstate_records.append({
            'receipt_id': receipt_id,
            'Molecule ID': row["Molecule ID"],
            'charge -4 (calc)': row["charge -4"],
            'charge -3 (calc)': row["charge -3"],
            'charge -2 (calc)': row["charge -2"],
            'charge -1 (calc)': row["charge -1"],
            'charge 0 (calc)': row["charge 0"],
            'charge 1 (calc)': row["charge 1"],
            'charge 2 (calc)': row["charge 2"],
            'charge 3 (calc)': row["charge 3"],
            'charge 4 (calc)': row["charge 4"],
            'charge -1 (exp)': exp_row["charge -1"].values[0],
            'charge 0 (exp)': exp_row["charge 0"].values[0],
            'charge 1 (exp)': exp_row["charge 1"].values[0],
            'charge 2 (exp)': exp_row["charge 2"].values[0],
        })
        
# Transform into Pandas DataFrame.
df_dom_ms_col = pd.DataFrame(data=dominant_microstate_records)
df_dom_ms_col
# Save dominant microstate collection as CSV
dom_ms_collection_filepath = "typeI_dominant_microstate_collection.csv"
df_dom_ms_col.to_csv(dom_ms_collection_filepath , index=False)

# Calculate the accuracy of dominant microstates
df_dom_ms_match = df_dom_ms_col.copy()
df_dom_ms_match  = df_dom_ms_match.drop(["charge -4 (calc)", "charge -3 (calc)",
                                      "charge -2 (calc)", "charge 3 (calc)", "charge 4 (calc)"], axis=1)

charges_of_exp_ms = [-1, 0, 1, 2]
for charge in charges_of_exp_ms:
    df_dom_ms_match["charge {} match".format(str(charge))] = np.NaN
df_dom_ms_match
# Populate the dataframe with 0 (False) or 1 (True) based on match of predicted microstates to experimental.
# Iterate through molecules
for i, row in enumerate(df_dom_ms_match.iterrows()):
    #print(i)
    #print(row)
    mol_ID = row[1]["Molecule ID"]
    # Iterate through charges
    for charge in charges_of_exp_ms:
        # Compare predicted and experimental microstate
        pred_ms = df_dom_ms_col.loc[df_dom_ms_col["Molecule ID"] == mol_ID]["charge {} (calc)".format(str(charge))].values[0]
        print("Pred:",pred_ms)
        exp_ms = df_dom_ms_col.loc[df_dom_ms_col["Molecule ID"] == mol_ID]["charge {} (exp)".format(str(charge))].values[0]
        print("Exp:", exp_ms)
        
        
        # if experimental microstate doesn't exist for that charge state match value should be NaN
        if isinstance(exp_ms, str) == False:
            if np.isnan(exp_ms):
                match = np.NaN
        # if experimental and predicted microstate IDs are the same, match value should be 1 (True)
        elif pred_ms == exp_ms:
            match = 1
        else:
            match = 0
            
        print("match: ", match)
        print()
        df_dom_ms_match.loc[df_dom_ms_col["Molecule ID"] == mol_ID , "charge {} match".format(str(charge))] = match
        #print(str(charge))
        #print(pred_row)
        #print(pred_row["charge {} (calc)".format(str(charge))])
    #["charge {} match".format{charge}]
    
df_dom_ms_match
# Calculate the accuracy of predicted microstate pairs<jupyter_output><empty_output><jupyter_text>2. Examples for Dominant Microstate Accuracy CalculationIf experimental microstate is not determined, match value will be NaN.If experiemental microstate exists but predicted microstate is missing or wrong, match value will be 0.If experimental and predicted microstate matches, match value will be 0<jupyter_code>pred_ms = np.NaN
exp_ms = "ms1"

# if experimental microstate doesn't exist for that charge state match value should be NaN
if isinstance(exp_ms, str) == False:
    if np.isnan(exp_ms):
        match = np.NaN
# if experimental and predicted microstate IDs are the same, match value should be 1 (True)
elif pred_ms == exp_ms:
    match = 1
else:
    match = 0
print(match)



pred_ms = "ms1"
exp_ms = np.NaN
# if experimental microstate doesn't exist for that charge state match value should be NaN
if isinstance(exp_ms, str) == False:
    if np.isnan(exp_ms):
        match = np.NaN
# if experimental and predicted microstate IDs are the same, match value should be 1 (True)
elif pred_ms == exp_ms:
    match = 1
else:
    match = 0
print(match)



pred_ms = "ms1"
exp_ms = "ms2"
# if experimental microstate doesn't exist for that charge state match value should be NaN
if isinstance(exp_ms, str) == False:
    if np.isnan(exp_ms):
        match = np.NaN
# if experimental and predicted microstate IDs are the same, match value should be 1 (True)
elif pred_ms == exp_ms:
    match = 1
else:
    match = 0
print(match)<jupyter_output>0
nan
0<jupyter_text>3. Analyze microstates of Jaguar and Epik<jupyter_code># Load prediction data
jaguartypei = SAMPL6DataProvider("typeI-raw-jaguar.csv", "typei", "Jaguar", bootstrap_options={"n_samples": 1})
epiktypei = SAMPL6DataProvider("typeI-raw-microscopic.csv", "typei", "Epik", bootstrap_options={"n_samples": 1})

# Load experimental data
exp = SAMPL6DataProvider("SAMPL6_experimental_pkas.csv", "exp", "Experiment", bootstrap_options={"n_samples": 3})
submission_IDs = ["Jaguar", "Epik"]
prediction_data_dict = {"Jaguar": jaguartypei, "Epik": epiktypei}
microstate_data = []

for receipt_id in submission_IDs:    
    for mol_ID in mol_IDs:
        predtypei = prediction_data_dict[receipt_id]
        pred_1mol = predtypei.load(mol_ID)
        microstate_IDs = pred_1mol.state_ids
        charges_of_microstates = pred_1mol.charges
        free_energies_of_microstates_pH0 = pred_1mol.free_energies[:,0]

#         print(mol_ID)
#         print("Microstates:\n", microstate_IDs)
#         print("Charges:", charges_of_microstates)
#         print("Free energies:\n",free_energies_of_microstates_pH0)
#         print()

        for i, microstate_ID in enumerate(microstate_IDs):

            microstate_data.append({
                'receipt_id': receipt_id,
                'Molecule ID': mol_ID,
                'Microstate ID': microstate_ID,
                'Charge': charges_of_microstates[i],
                '$\Delta$G (pH=0)': free_energies_of_microstates_pH0[i]
            })

df_microstate_data = pd.DataFrame(data=microstate_data)
df_microstate_data
#pred_mol_IDs = set(df_microstate_data["Molecule ID"].values)
#pred_mol_IDs

df_jaguar = df_microstate_data[df_microstate_data["receipt_id"] == "Jaguar"]
print(df_jaguar[df_jaguar["Molecule ID"] == "SM07"])

df_epik = df_microstate_data[df_microstate_data["receipt_id"] == "Epik"]
print(df_epik[df_epik["Molecule ID"] == "SM07"])

# Create an empty dataframe to store dominant predicted microstate ID of each charge 

dominant_microstate_records = []

for receipt_id in submission_IDs:
    for mol_ID in pred_mol_IDs:     
        # Append to a list to later save as a CSV
        dominant_microstate_records.append({
            'receipt_id': receipt_id,
            'Molecule ID': mol_ID,
            'charge -4': np.NaN,
            'charge -3': np.NaN,
            'charge -2': np.NaN,
            'charge -1': np.NaN,
            'charge 0': np.NaN,
            'charge 1': np.NaN,
            'charge 2': np.NaN,
            'charge 3': np.NaN,
            'charge 4': np.NaN
        })
        
# Transform into Pandas DataFrame.
df_dominant_ms = pd.DataFrame(data=dominant_microstate_records)
df_dominant_ms
# Populate data frame with predicted dominant microstates. 
# If prediction of a charge state is missing it will be recorded as NaN

for receipt_id in submission_IDs:
    
    # Get set fo pred_mol_IDs from the predicted microstates table
    df_1submission = df_microstate_data[df_microstate_data["receipt_id"] == receipt_id]
    #print("df_1submission:\n", df_1submission.head())
    
    #df_SM07 = df_1submission[df_1submission["Molecule ID"] == "SM07"]
    #df_SM07_charge1 = df_SM07[df_SM07["Charge"] == 1]
    #print("df_1submission for SM07 charge 1:\n", df_SM07_charge1)
    #idxmin = df_SM07_charge1['$\Delta$G (pH=0)'].idxmin()
    #print("idxmin:", idxmin)
    #SM07_charge1_dom_ms = df_SM07_charge1.loc[idxmin, "Microstate ID"]
    #print("SM07_charge1_dom_ms:", SM07_charge1_dom_ms)
    
    pred_mol_IDs = set(df_1submission["Molecule ID"].values)
    
    # Record dominant microstate series for each molecule
    for mol_ID in pred_mol_IDs:     
        df_1mol = df_1submission[df_1submission["Molecule ID"] == mol_ID]
        charges = set(df_1mol["Charge"].values)

        for i, charge in enumerate(charges):
            df_1mol_1charge = df_1mol[df_1mol["Charge"] == charge]
            #dominant_microstate = df_1mol_1charge.loc[df_1mol_1charge['$\Delta$G (pH=0)'].idxmin()]["Microstate ID"]
            dominant_microstate = df_1mol_1charge.loc[df_1mol_1charge['$\Delta$G (pH=0)'].idxmin(),"Microstate ID"]
            df_dominant_ms.loc[(df_dominant_ms["Molecule ID"] == mol_ID) & (df_dominant_ms["receipt_id"] == receipt_id), "charge {}".format(charge)] = dominant_microstate

            #charge_and_dom_ms_dict.update({charge : dominant_microstate})

# Save as CSV
df_dominant_ms.to_csv("typeI_submission_dominant_microstate_collection_only_predictions.csv", index=False)

df_dominant_ms
# Read experimental dominant microstates data
exp_microstates = pd.read_csv("experimental_microstates_with_charge.csv")

# Organize experimental microstate in dominant microstate collection format

# Create empty dataframe
df_exp_dom_ms = pd.DataFrame(columns=[ "Molecule ID", "charge -4", "charge -3","charge -2", "charge -1","charge 0",
                                       "charge 1","charge 2", "charge 3", "charge 4"])
charges = np.arange(-4,4,1)
exp_mol_IDs = set(exp_microstates["Molecule ID"].values)

for i, mol_ID in enumerate(exp_mol_IDs): 
    df_exp_dom_ms.loc[i] = [mol_ID, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN, np.NaN ]

# Make molecule ID the index
df_exp_dom_ms = df_exp_dom_ms.set_index("Molecule ID", drop = False)

# Populate expermental dominant microstate collection
for i, mol_ID in enumerate(exp_mol_IDs): 
    for charge in charges:
        #print("mol {}, charge {}".format(mol_ID, charge))
        exp_ms_1mol = exp_microstates[exp_microstates["Molecule ID"] == mol_ID]
        #print("exp_ms_1mol:\n",exp_ms_1mol)
        
        # Check if charge exist in "Charge of A" column
        exp_charge_of_A = exp_ms_1mol["Charge of A"].values
        exp_charge_of_HA = exp_ms_1mol["Charge of HA"].values
        #print("exp_charge_of_A:", exp_charge_of_A)
        #print("exp_charge_of_HA:", exp_charge_of_HA)
        
        if charge in exp_charge_of_A:
            #print("Charges match A.")
            #exp_dom_ms_ID =  exp_ms_1mol["Microstate ID of A"].values
            exp_dom_ms_ID = exp_ms_1mol[exp_ms_1mol["Charge of A"] == charge]["Microstate ID of A"].values[0]
        elif charge in exp_charge_of_HA:
            #print("Charges match HA.")
            #exp_dom_ms_ID =  exp_ms_1mol["Microstate ID of HA"].values
            exp_dom_ms_ID = exp_ms_1mol[exp_ms_1mol["Charge of HA"] == charge]["Microstate ID of HA"].values[0]
        else:
            #print("Different charge")
            exp_dom_ms_ID = np.NaN
        
        #print("exp_dom_ms_ID:", exp_dom_ms_ID)
        
        df_exp_dom_ms.loc[mol_ID, "charge {}".format(charge)] = exp_dom_ms_ID

df_exp_dom_ms.to_csv("typeI_submission_dominant_microstate_collection_only_experiments.csv", index=False)
df_exp_dom_ms
# Create a collection file with both experimental and predicted microstate IDs

# Read experimental data collection
df_exp_dom_ms = pd.read_csv("typeI_submission_dominant_microstate_collection_only_experiments.csv")
                           
# Read predicted microstates collections
df_pred_dom_ms = pd.read_csv("typeI_submission_dominant_microstate_collection_only_predictions.csv")

# Record data in a list to later convert to collection dataframe
dominant_microstate_records = []

receipt_id_list = set(df_pred_dom_ms["receipt_id"].values)

for receipt_id in receipt_id_list:
    df_1method = df_pred_dom_ms[df_pred_dom_ms.receipt_id == receipt_id]
    for i, row in df_1method.iterrows():
        mol_ID = row["Molecule ID"]
        
        # Find the matching row from experimental dominant microstate dataframe using the mol_ID
        exp_row = df_exp_dom_ms[df_exp_dom_ms["Molecule ID"] == mol_ID]
        #print(exp_row)
        mol_ID_exp = exp_row["Molecule ID"].values[0]
        
        
        # Append to a list to later save as a CSV
        dominant_microstate_records.append({
            'receipt_id': receipt_id,
            'Molecule ID': row["Molecule ID"],
            'charge -4 (calc)': row["charge -4"],
            'charge -3 (calc)': row["charge -3"],
            'charge -2 (calc)': row["charge -2"],
            'charge -1 (calc)': row["charge -1"],
            'charge 0 (calc)': row["charge 0"],
            'charge 1 (calc)': row["charge 1"],
            'charge 2 (calc)': row["charge 2"],
            'charge 3 (calc)': row["charge 3"],
            'charge 4 (calc)': row["charge 4"],
            'charge -1 (exp)': exp_row["charge -1"].values[0],
            'charge 0 (exp)': exp_row["charge 0"].values[0],
            'charge 1 (exp)': exp_row["charge 1"].values[0],
            'charge 2 (exp)': exp_row["charge 2"].values[0],
        })
        
# Transform into Pandas DataFrame.
df_dom_ms_col = pd.DataFrame(data=dominant_microstate_records)
df_dom_ms_col
# Calculate the accuracy of dominant microstates
# Experimental data doesn't include any microstates for charges -4, -3, -1, +3, and +4. 
# So these won't be analyzed for accuracy.
df_dom_ms_match = df_dom_ms_col.copy()
df_dom_ms_match  = df_dom_ms_match.drop(["charge -4 (calc)", "charge -3 (calc)",
                                      "charge -2 (calc)", "charge 3 (calc)", "charge 4 (calc)"], axis=1)
charges_of_exp_ms = [-1, 0, 1, 2]

# Create empty columns to record match accuracy
for charge in charges_of_exp_ms:
    df_dom_ms_match["charge {} match".format(str(charge))] = np.NaN

# Populate the dataframe with 0 (False) or 1 (True) based on match of predicted microstates to experimental.
# Iterate through submissions
receipt_id_list = set(df_dom_ms_match["receipt_id"].values)

for receipt_id in receipt_id_list:
    df_1submission = df_dom_ms_match[df_dom_ms_match["receipt_id"] == receipt_id]

    # Iterate through molecules
    for i, row in enumerate(df_1submission.iterrows()):
        mol_ID = row[1]["Molecule ID"]
        # Iterate through charges

        for charge in charges_of_exp_ms:
            # Compare predicted and experimental microstate
            pred_ms = df_1submission.loc[df_1submission["Molecule ID"] == mol_ID]["charge {} (calc)".format(str(charge))].values[0]
            #print("Pred:",pred_ms)
            exp_ms = df_1submission.loc[df_1submission["Molecule ID"] == mol_ID]["charge {} (exp)".format(str(charge))].values[0]
            #print("Exp:", exp_ms)

            # if experimental microstate doesn't exist for that charge state match value should be NaN
            if isinstance(exp_ms, str) == False:
                if np.isnan(exp_ms):
                    match = np.NaN
            # if experimental and predicted microstate IDs are the same, match value should be 1 (True)
            elif pred_ms == exp_ms:
                match = 1
            else:
                match = 0
            #print("match: ", match)
            #print()
            
            # Record the match value to the dominant microstates collection dataframe
            df_dom_ms_match.loc[(df_dom_ms_col["Molecule ID"] == mol_ID) & (df_dom_ms_col["receipt_id"] == receipt_id) , "charge {} match".format(str(charge))] = match
            #print(str(charge))
            #print(pred_row)
            #print(pred_row["charge {} (calc)".format(str(charge))])

# Save dominant microstate collection as CSV
dom_ms_collection_filepath = "typeI_dominant_microstate_collection.csv"
df_dom_ms_col.to_csv(dom_ms_collection_filepath , index=False)

df_dom_ms_match<jupyter_output><empty_output><jupyter_text>3.1 Calculate overall dominant microstate match accuracy for each method<jupyter_code># Calculate overall dominant microstate match accuracy for each method

# Save data in a list of dictionaries to convert to dataframe later
dom_ms_stats = []

for receipt_id in submission_IDs:
    # Take subset of dominant microstate collection based on submission ID
    df_1submission = df_dom_ms_match[df_dom_ms_match["receipt_id"] == receipt_id]
    
    match_values = []
    for charge in charges_of_exp_ms:
        match_values_of_1charge = df_1submission["charge {} match".format(charge)].values
        #print(match_values_of_1charge)
        for value in match_values_of_1charge:
            match_values.append(value)
    match_values_array = np.array(match_values)       
    print(match_values_array)
    
    # Calculate overall dominant microstate match accuracy ignoring NaNs
    accuracy = np.nanmean(match_values_array)
    print("Accuracy:", accuracy)
    
    dom_ms_stats.append({
            'receipt_id': receipt_id,
            'Dominant Microstate Accuracy': accuracy
    })

df_dom_ms_stats = pd.DataFrame(data=dom_ms_stats)
df_dom_ms_stats 
import matplotlib
import matplotlib.pyplot as plt
import seaborn as sns
%matplotlib inline
plt.rcParams["figure.figsize"] = [4,4]

def barplot(df, x_label, y_label, title):
    """Creates bar plot of a given dataframe.

    Args:
        df: Pandas Dataframe that should have columns with columnnames specified in other arguments.
        x_label: str, column name of x axis categories
        y_label: str, column name of y axis values
        title: str, the title of the plot

    """
    # Plot style
    plt.close()
    plt.style.use(["seaborn-talk", "seaborn-whitegrid"])
    plt.rcParams['axes.labelsize'] = 18
    plt.rcParams['xtick.labelsize'] = 14
    plt.rcParams['ytick.labelsize'] = 16
    plt.ylim(0,1)
    #plt.tight_layout()

    # Plot
    data = df
    x = range(len(data[y_label]))
    y = data[y_label]
    plt.bar(x, y)
    plt.xticks(x, data[x_label], rotation=90)
    plt.xlabel(x_label)
    plt.ylabel(y_label)
    if len(title) > 70:
        plt.title(title[:70]+"...")
    else:
        plt.title(title)
    plt.tight_layout()
    

def barplot_with_2groups(df, x_label, y_label1, y_label2):
    """Creates bar plot of a given dataframe with 2 datasets.

    Args:
        df: Pandas Dataframe that should have columns with columnnames specified in other arguments.
        x_label: str, column name of x axis categories
        y_label1: str, column name of y axis values 
        y_lower2: str, column name of y axis values
      
    """

    # Color
    current_palette = sns.color_palette()
    #current_palette = sns.color_palette("GnBu_d")
    error_color = sns.color_palette("GnBu_d")[0]

    # Plot style
    plt.close()
    plt.style.use(["seaborn-talk", "seaborn-whitegrid"])
    plt.rcParams['axes.labelsize'] = 18
    plt.rcParams['xtick.labelsize'] = 14
    plt.rcParams['ytick.labelsize'] = 16
    plt.tight_layout()
    bar_width = 0.45

    # Plot 1st group of data
    data = df  # Pandas DataFrame    
    x = range(len(data[x_label]))
    y = data[y_label1]
    plt.bar(x, y, label = y_label1, width=bar_width, color=current_palette[0])
    plt.xticks(x, data[x_label], rotation=90)
    
    # Plot 2nd group of data
    data = df  # Pandas DataFrame
    index = np.arange(data.shape[0])

    x = range(len(data[x_label]))
    y = data[y_label2]
    #plt.bar(x, y)
    plt.bar(index + bar_width, y, label = y_label2, width=bar_width, color=sns.color_palette("BuGn_r")[3])
    plt.xticks(index + bar_width/2, data[x_label], rotation=90)
    
    plt.xlabel(x_label)
    plt.ylabel("Dominant Microstate Accuracy")
    plt.legend( bbox_to_anchor=(1.05, 0, 0.5, 1))
df_dom_ms_stats_sorted = df_dom_ms_stats.sort_values(by="Dominant Microstate Accuracy", ascending=False)

barplot(df = df_dom_ms_stats_sorted, 
        x_label = "receipt_id", 
        y_label = "Dominant Microstate Accuracy", 
        title = "Accuracy of dominant microstate predictions in Type I submissions" )
plt.savefig("typeI_dominant_microstate_accuracy.pdf")
# Calculate dominant microstate match accuracy for each method for neutral and +1 charge state separately

charges_to_be_calculated = [0, 1]

for charge in charges_to_be_calculated:
    df_dom_ms_stats["Dominant Microstate Accuracy of Charge {}".format(str(charge))] = np.NaN

    for receipt_id in submission_IDs:
        # Take subset of dominant microstate collection based on submission ID
        df_1submission = df_dom_ms_match[df_dom_ms_match["receipt_id"] == receipt_id]

        match_values_of_1charge = df_1submission["charge {} match".format(charge)].values
        match_values_array = np.array(match_values_of_1charge)       

        # Calculate overall dominant microstate match accuracy ignoring NaNs
        accuracy = np.nanmean(match_values_array)
        
        # Record accuracy value to the dataframe
        df_dom_ms_stats.loc[df_dom_ms_stats["receipt_id"] == receipt_id, "Dominant Microstate Accuracy of Charge {}".format(str(charge))] = accuracy
        
df_dom_ms_stats
barplot_with_2groups(df = df_dom_ms_stats, 
                     x_label = "receipt_id", 
                     y_label1 = "Dominant Microstate Accuracy of Charge 0", 
                     y_label2 = "Dominant Microstate Accuracy of Charge 1")
plt.savefig("typeI_dominant_microstate_accuracy_separate_charges.pdf")<jupyter_output><empty_output><jupyter_text>3.2 Calculate overall dominant microstate match accuracy for each molecule<jupyter_code># Plot dominant microstate error vs Molecule ID

# Save data in a list of dictionaries to convert to dataframe later
dom_ms_stats = []

mol_IDs = set(df_dom_ms_match["Molecule ID"].values)

for mol_ID in mol_IDs:
    # Take subset of dominant microstate collection based on molecule ID
    df_1molecule = df_dom_ms_match[df_dom_ms_match["Molecule ID"] == mol_ID]
    #print(mol_ID)
    #print("df_1molecule:\n", df_1molecule)
    
    match_values = []
    for charge in charges_of_exp_ms:
        match_values_of_1charge = df_1molecule["charge {} match".format(charge)].values
        #print(match_values_of_1charge)
        for value in match_values_of_1charge:
            match_values.append(value)
    match_values_array = np.array(match_values)       
    #print(match_values_array)

    # Calculate overall dominant microstate match accuracy ignoring NaNs
    accuracy = np.nanmean(match_values_array)
    #print("Accuracy:", accuracy)

    dom_ms_stats.append({
            'Molecule ID': mol_ID,
            'Dominant Microstate Accuracy': accuracy
    })

df_dom_ms_stats_mol_wise = pd.DataFrame(data=dom_ms_stats)
df_dom_ms_stats_mol_wise
df_dom_ms_stats_mol_wise_sorted = df_dom_ms_stats_mol_wise.sort_values(by="Dominant Microstate Accuracy", ascending=False)

barplot(df = df_dom_ms_stats_mol_wise_sorted, 
        x_label = "Molecule ID", 
        y_label = "Dominant Microstate Accuracy", 
        title = "Accuracy of dominant microstate predictions in Type I submissions" )
plt.savefig("typeI_dominant_microstate_accuracy_for_each_molecule.pdf")
# Calculate dominant microstate match accuracy for each molecule for neutral and +1 charge state separately

charges_to_be_calculated = [0, 1]

for charge in charges_to_be_calculated:
    df_dom_ms_stats_mol_wise["Dominant Microstate Accuracy of Charge {}".format(str(charge))] = np.NaN

    for mol_ID in mol_IDs:
        # Take subset of dominant microstate collection based on submission ID
        df_1molecule = df_dom_ms_match[df_dom_ms_match["Molecule ID"] == mol_ID]

        match_values_of_1charge = df_1molecule["charge {} match".format(charge)].values
        match_values_array = np.array(match_values_of_1charge)       

        # Calculate overall dominant microstate match accuracy ignoring NaNs
        accuracy = np.nanmean(match_values_array)
        
        # Record accuracy value to the dataframe
        df_dom_ms_stats_mol_wise.loc[df_dom_ms_stats_mol_wise["Molecule ID"] == mol_ID, "Dominant Microstate Accuracy of Charge {}".format(str(charge))] = accuracy
        
df_dom_ms_stats_mol_wise
barplot_with_2groups(df = df_dom_ms_stats_mol_wise, 
                     x_label = "Molecule ID", 
                     y_label1 = "Dominant Microstate Accuracy of Charge 0", 
                     y_label2 = "Dominant Microstate Accuracy of Charge 1")
plt.savefig("typeI_dominant_microstate_accuracy_for_each_molecule_separate_charges.pdf")<jupyter_output><empty_output><jupyter_text>3.3 Calculate overall dominant microstate match accuracy for each method with bootstrapped confidence intervals<jupyter_code>def nanmean(data):
    #x, y = data.T
    return np.nanmean(data)

def compute_bootstrap_statistics(samples, stats_funcs, percentile=0.95, n_bootstrap_samples=10000): #10000
    """Compute bootstrap confidence interval for the given statistics functions."""
    # Handle case where only a single function is passed.
    try:
        len(stats_funcs)
    except TypeError:
        stats_funcs = [stats_funcs]

    # Compute mean statistics.
    statistics = [stats_func(samples) for stats_func in stats_funcs]

    # Generate bootstrap statistics.
    bootstrap_samples_statistics = np.zeros((len(statistics), n_bootstrap_samples))
    for bootstrap_sample_idx in range(n_bootstrap_samples):
        samples_indices = np.random.randint(low=0, high=len(samples), size=len(samples))
        for stats_func_idx, stats_func in enumerate(stats_funcs):
            bootstrap_samples_statistics[stats_func_idx][bootstrap_sample_idx] = stats_func(samples[samples_indices])

    # Compute confidence intervals.
    percentile_index = int(np.floor(n_bootstrap_samples * (1 - percentile) / 2)) - 1
    bootstrap_statistics = []
    for stats_func_idx, samples_statistics in enumerate(bootstrap_samples_statistics):
        samples_statistics.sort()
        stat_lower_percentile = samples_statistics[percentile_index]
        stat_higher_percentile = samples_statistics[-percentile_index+1]
        confidence_interval = (stat_lower_percentile, stat_higher_percentile)
        bootstrap_statistics.append([statistics[stats_func_idx], confidence_interval, samples_statistics])

    return bootstrap_statistics

#stats_funcs = collections.OrderedDict([('Accuracy', nanmean)])
#stats_funcs = nanmean

#bootstrap_statistics = compute_bootstrap_statistics(data.as_matrix(), stats_funcs,
#                                                                n_bootstrap_samples=10000)

# Calculate overall dominant microstate match accuracy for each method

# Save data in a list of dictionaries to convert to dataframe later
dom_ms_stats = []

for receipt_id in submission_IDs:
    # Take subset of dominant microstate collection based on submission ID
    df_1submission = df_dom_ms_match[df_dom_ms_match["receipt_id"] == receipt_id]
    
    match_values = []
    for charge in charges_of_exp_ms:
        match_values_of_1charge = df_1submission["charge {} match".format(charge)].values
        #print(match_values_of_1charge)
        for value in match_values_of_1charge:
            match_values.append(value)
    match_values_array = np.array(match_values)       
    print(match_values_array)
    
    # Calculate overall dominant microstate match accuracy ignoring NaNs
    bootstrap_statistics = compute_bootstrap_statistics(match_values_array, stats_funcs = nanmean,
                                                                n_bootstrap_samples=1000)
    #print("bootstrap_statistics:\n", bootstrap_statistics)
    accuracy = bootstrap_statistics[0][0]
    accuracy_lower_bound, accuracy_upper_bound = bootstrap_statistics[0][1]
    
    #accuracy = np.nanmean(match_values_array)
    print("Accuracy:", accuracy)
    print("CI: ", accuracy_lower_bound, accuracy_upper_bound)
    
    dom_ms_stats.append({
            'receipt_id': receipt_id,
            'Accuracy': accuracy,
            'Accuracy lower bound': accuracy_lower_bound,
            'Accuracy upper bound': accuracy_upper_bound
    })

df_dom_ms_stats = pd.DataFrame(data=dom_ms_stats)
df_dom_ms_stats<jupyter_output>[ 1. nan nan nan nan nan nan nan  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.
  1.  1.  1.  1.  1.  1. nan nan nan nan nan  1. nan nan]
Accuracy: 1.0
CI:  1.0 1.0
[ 1. nan nan nan nan nan nan nan  1.  1.  1.  1.  1.  1.  1.  1.  1.  0.
  1.  1.  1.  1.  0.  1. nan nan nan nan nan  1. nan nan]
Accuracy: 0.8888888888888888
CI:  0.7222222222222222 1.0
{'fraction_non_alphanumeric': 0.10507969527450595, 'fraction_numerical': 0.015433053158294211, 'mean_word_length': 3.4545454545454546, 'pattern_counts': {'":': 7, '<': 23, '<?xml version=': 0, '>': 24, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 2, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '344339', 'license': 'MIT', 'n_tokens_mistral': 12309, 'n_tokens_neox': 11555, 'n_words': 2731, 'path': 'analysis_of_pKa_predictions/20191003_typeI_microstate_deltaG_calculation/calc_relative_free_energy_at_pH_0.ipynb', 'repo_name': 'choderalab/sampl6-logD-compound-selection', 'size': '35573'}