Mercurial > repos > jay > gaiac_violin_plot
comparison gaiac_bias_calculation/gaiac_bias_calculation.py @ 0:abfc2c9779d6 draft
planemo upload for repository https://github.com/jaidevjoshi83/gaiac.git commit c29a769ed165f313a6410925be24f776652a9663-dirty
| author | jay |
|---|---|
| date | Thu, 15 May 2025 14:47:54 +0000 |
| parents | |
| children |
comparison
equal
deleted
inserted
replaced
| -1:000000000000 | 0:abfc2c9779d6 |
|---|---|
| 1 import scipy | |
| 2 #from scipy.stats.distributions import chi2,t | |
| 3 from scipy.stats import t | |
| 4 from scipy import stats | |
| 5 import pandas as pd | |
| 6 import numpy as np | |
| 7 import argparse | |
| 8 | |
| 9 def Bias_abs(infile, clm1, clm2, Out): | |
| 10 | |
| 11 df = pd.read_csv(infile, sep="\t") | |
| 12 lOc = int(df.shape[0]) | |
| 13 cl = df.columns.tolist() | |
| 14 #making a column of absolute di as per us epa | |
| 15 df['PM10_OPC_DIFF']=abs(((df[cl[int(clm1)-1]])-(df[cl[int(clm2)-1]]))/(df[cl[int(clm2)-1]]))*100 | |
| 16 #square of di | |
| 17 df['PM10_OPC_DIFFs']=(df['PM10_OPC_DIFF'])*(df['PM10_OPC_DIFF']) | |
| 18 #summation of the columns | |
| 19 d10=df.PM10_OPC_DIFF.sum() | |
| 20 d10_2=df.PM10_OPC_DIFFs.sum() | |
| 21 #AB and AS calculations | |
| 22 AB= d10/lOc | |
| 23 AS=(np.sqrt((lOc*(d10_2)-(d10)**2)/(2*lOc*(lOc-1)))) | |
| 24 | |
| 25 #T distribution calculation | |
| 26 T=stats.t.ppf(1-0.05, lOc) | |
| 27 #Absolute bias calculation | |
| 28 Bias_abs = abs(AB + T*AS/(np.sqrt(lOc))) | |
| 29 | |
| 30 #di column with sign | |
| 31 df['Di']=(((df[cl[int(clm1)-1]])-(df[cl[int(clm2)-1]]))/(df[cl[int(clm2)-1]]))*100 | |
| 32 #quantiles of di without sign | |
| 33 q1=df['Di'].quantile([0.25]) | |
| 34 q2=df['Di'].quantile([0.75]) | |
| 35 #assigning sign to absolute bias based on q1 and q2 values | |
| 36 | |
| 37 if (q1[0.25] < 0) & (q2[0.75] < 0): | |
| 38 Bias_abs=Bias_abs*-1 | |
| 39 elif (q1[0.25] > 0) & (q2[0.75] > 0): | |
| 40 Bias_abs=Bias_abs | |
| 41 elif (q1[0.25] >0) & (q2[0.75]<0): | |
| 42 Bias_abs=(u"\u00B1"+str(Bias_abs)) | |
| 43 elif (q1[0.25] <0) & (q2[0.75]>0): | |
| 44 Bias_abs=(u"\u00B1"+str(Bias_abs)) | |
| 45 | |
| 46 #output file | |
| 47 df1 = pd.DataFrame([Bias_abs], columns=['Percent Bias']) | |
| 48 df1.round(4).to_csv(Out, sep="\t") | |
| 49 | |
| 50 | |
| 51 if __name__=="__main__": | |
| 52 | |
| 53 parser = argparse.ArgumentParser() | |
| 54 | |
| 55 parser.add_argument("-I", "--infile", required=True, default=None, help="Input file") | |
| 56 parser.add_argument("-c1", "--column_1", required=True, default=None, help="First column") | |
| 57 parser.add_argument("-c2", "--column_2", required=True, default=None, help="Second column") | |
| 58 parser.add_argument("-o", "--output", required=True, default=None, help="OutFile") | |
| 59 args = parser.parse_args() | |
| 60 | |
| 61 Bias_abs(args.infile, args.column_1, args.column_2, args.output) |
