This is my first script that I've written. As a result, I'm sure there are extra lines that are unneeded, or maybe better more concise ways of doing things than I have done here. I have tried to add explanations to each section so that it's a bit easier to understand what's going on, and why I'm doing what I'm doing.
Just a bit of backstory to explain what this script is designed to do. This script is designed to take in the output files from a specific program that is used to analyze we'll call "Titration Data", and to extract the data from these files and do some modifications to them. Once they've been modified, it then tries to fit the data to a model, and output the graphs, fit value itself, standard deviation, goodness of fit, and the output of the function (so someone can graph it in another program if they want). There are only 3 input files that are needed: The data files themselves, the concentration file, and the dilution file. The dimensions for all of these need to match.
That said, below is my script, any advice/feedback is greatly appreciated! If there is anything unclear, please let me know and I can clear it up. Thank you!
#I plan on adding other scripts later to analyze different data sets. Currently I have written the script for one. But this is designed to be a general use program, so I want to get it out there now, but add features to it later. Hence all the initial inputs. A=input('Slow or Fast Exchange? (Type in "Slow" or "Fast")\n') if A == ('Slow').lower(): B=input('Will you be using both the free and bound state peaks? (Type in "Yes" or "No")\n') if B==('Yes').lower(): print ('Slow Exchange analysis using both states is not developed yet') else: C=input('Will you be using the free or bound state peaks? (Type in "Free" or "Bound")\n') if C==('Free').lower(): print ('Loading up Slow Exchange Script') import re import os import pandas as pd import numpy as np import matplotlib.pyplot as plt from scipy.optimize import curve_fit #This is where the user inputs where the files they will be calling are stored. To make this more general user friendly, I automatically fill the drive and user portion main_path=(r'C:\Users\Sams PC') next_path=input('Please indicate the folder the files are. E.g. /Desktop/Peaklists \n') full_path=main_path+next_path os.chdir(full_path) #The data that will be analyzed is split into 3 sections. One is the experimental data itself (Titration_files), the settings they used to obtain the experimental data (Concentrations), and a normalization factor (Dilution) Titration_files =  Data_Table =  Peak_Height =  Titration_file_input = '' Peak_Names= Dilutions=input('Enter name of Dilution file.\n') Concentrations=input('Enter name of Concentration file.\n') # This is designed so the user can input in as many files as they desire. And if they misspell, or accidently add an extra file or forget a file, the script will still continue. print('Enter name of peaklist files. When finished, type done and enter to stop.') while True: Titration_file_input = input() if Titration_file_input.lower() not in ['done']: Titration_files.append(Titration_file_input) else: break for Data in Titration_files: try: Titration_Datatable = pd.read_csv(Data, sep='\s+', header=None) Titration_Datatable.columns=['Column_1','Column_2','Column_3', 'Column_4', 'Column_5'] Data_Table.append(Titration_Datatable) except: print('File' + ' ' + Data + ' ' + 'not found') #The input of Titration_Data comes from a specific program that has the required data stored in Column 4. I am creating a second table from Column 1 to use when saving the plots and savefiles. for Titration_Datatable in Data_Table: Peak_Height.append(Titration_Datatable.loc[:,'Column_4']) Peak_Names.append(Titration_Datatable.loc[:,'Column_1'].drop(,axis=0).drop(,axis=0)) #The first row in the datatable automatically is labeled Data, so that is removed. The 2nd row also has miscellanious info. The .astype is because I had to define them as integers, otherwise I would get errors at the division step below. concatenated_Titration_Datatable = pd.concat(Peak_Height, axis=1).drop(,axis=0).drop(,axis=0).astype(int) Dilutions=pd.read_csv('Dilution.txt', sep='\s+', header=None) #The Data may sometimes have negative values, which can be considered the same as zero in this case. Combined=concatenated_Titration_Datatable.clip(lower=0) #The Data now needs to be normalized to account for Dilution (this is introduced in the experimental setup), the values are defined by a txt file the user uploads. Normalized=(Combined/Dilutions.values) #It appears easier to do data modifications with a numpy matrix rather than with pandas, so I converted this to numpy. I also had to do it for the function below as well. M=pd.DataFrame.to_numpy(Normalized) #The below function is part of Data processing. It's an equation for analyzing the data. J=(M[:, :1]-M)/((M[:, :1]-M)+ M) Titration_Data=J Concentrations=np.loadtxt('Concentrations.txt') Protein=Concentrations[:,0] Ligand=Concentrations[:,1] Input_Data=[Protein,Ligand] x=Input_Data A=Input_Data B=Input_Data+Input_Data C=Input_Data #To be able to save the Peak_Names list above as a png for the graphs below, I had to convert to a nupy array first. But this would cause errors since its a matrix, thus I removed all the repeats, giving me a 1D array. Peak_Names_ar=np.array(Peak_Names) Peak_Names_array=np.unique(Peak_Names_ar) #I have setup the values A,B, and C because it makes the below equation cleaner and easier. The only reason x is defined is because I couldn't run the function without it, although it doesn't appear to effect the output. def fun(x, kd): return np.array((B+kd-np.sqrt(((B+kd)**2)-4*A*C))/(2*A)) #For each iteration of the function, I want to save the output of the fit (kD), the standard deviation of that fit, and R2 (goodness of fit). I'm also saving each plot, and using the above Peak_Names file to do so. kD= r2= standard_deviation= output_for_graphing= for values,i in zip(Titration_Data,Peak_Names_array): Intensity=[values] Intensity_Array=np.array(Intensity) y=Intensity_Array.flatten() popt, pcov = curve_fit(fun, x, y) kD.append(popt) fun_data=fun(x,*popt) output_for_graphing.append(fun_data) residuals=y-fun(x, popt) ss_res=np.sum(residuals**2) ss_tot=np.sum((y-np.mean(y))**2) r_squared=1-(ss_res/ss_tot) r2.append(r_squared) std = np.sqrt(np.diag(pcov)) standard_deviation.append(std) plt.plot(x, y, label='data') plt.plot(x, fun(x, *popt), label='fitted') plt.xlabel('Ligand Concentration') plt.ylabel('Intensity') plt.title([i]) plt.grid() plt.legend() files_to_save=str([i])+'.png' plt.savefig(files_to_save) plt.show() #I'm saving the output of the function as well, in case someone wants to graph it for themselves and doesn't like matplotlibs display. np.savetxt('Output for graphing.txt',output_for_graphing) #I need to flatten the arrays to be able to stack them (I get an error otherwise) kD_array=np.array(kD).flatten() r2_array=np.array(r2).flatten() standard_deviation_array=np.array(standard_deviation).flatten() Dissociation_Constant=np.stack((kD_array,standard_deviation_array,r2_array),axis=-1) # I realized after the fact that I wanted to add an extra column to the above matrix. But since numpy doesn't seem to have an insert command, and I didn't want to break the matrix I already had, I thought it might be easier to just change it to a pandas datatable. Dissociaton_Constant_Table=pd.DataFrame(Dissociation_Constant) Dissociaton_Constant_Table.columns=['kD', 'Standard Deviation', 'R2'] Dissociaton_Constant_Table.insert(0,'Files',Peak_Names_array) Dissociaton_Constant_Table.to_csv('Dissociation_Constant.txt', sep='\t', index=False) else: print ('Slow exchange analysis for the bound peak is not developed yet')
Input Files (incase anyone wanted to run it themselves):
#Concentrations 0.6 0.06 0.596421471 0.119284294 0.5859375 0.29296875 #Dilution 1.03E+05 1.03E+05 9.50E+04 #Data file 1 Assignment w1 w2 Data Height Dil_XN-HN 122.352 7.682 102832 1XH-HN 111.783 9.206 50034 2XN-HN 104.093 9.040 38230 3XN-HN 107.595 7.443 34503 4XN-HN 106.481 8.252 34219 5XN-HN 108.303 8.396 37442 6XN-HN 112.540 9.481 32144 7XN-HN 132.030 9.051 28053 8XN-HN 128.895 8.000 30395 #Data file 2 Assignment w1 w2 Data Height Dil_XN-HN 122.352 7.682 102832 1XH-HN 111.771 9.206 50034 2XN-HN 104.093 9.040 38230 3XN-HN 107.595 7.443 34503 4XN-HN 106.481 8.252 34219 5XN-HN 108.303 8.396 37442 6XN-HN 112.540 9.481 32144 7XN-HN 132.030 9.051 28053 8XN-HN 128.895 8.000 30395 #Data file 3: Assignment w1 w2 Data Height Dil_XN-HN 122.348 7.679 95020 1XH-HN 111.738 9.202 32355 2XN-HN 104.118 9.037 31054 3XN-HN 107.619 7.438 22208 4XN-HN 106.516 8.252 22975 5XN-HN 108.300 8.397 25319 6XN-HN 112.506 9.484 21293 7XN-HN 132.037 9.047 19972 8XN-HN 128.906 7.992 22899 ```