14-03-2021 Guido Cattani
from pathlib import Path
import pandas as pd
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
import matplotlib.lines as mlines
import matplotlib.patches as mpatches
from scipy.stats import linregress as linregress
# read difference sensation levels BC - AC simulated with a 65 dB ISDS input signal at BCD side (90 degree angle)
def read_diff_sl():
f_in = '/media/guido/LACIE/Cingle_Guido/Master/Diff_SL.xlsx'
p_in = Path(f_in)
df = pd.read_excel(p_in, engine='openpyxl', sheet_name='Diff_SL_65dB_90deg', header=0, nrows=85)
df = df.drop(['Unnamed: 0'], axis=1)
df = df.drop(columns='Study_ID')
df = df.fillna(pd.NA)
return df
# read measured thresholds
def read_AC_thresholds():
f_in = '/media/guido/LACIE/Cingle_Guido/Master/AC_thresholds.xlsx'
p_in = Path(f_in)
df = pd.read_excel(p_in, engine='openpyxl', sheet_name=0, header=0, nrows=85)
df = df.drop(['Unnamed: 0'], axis=1)
type_dict = {'Study_ID': 'str', 'Device': 'str',
'AC_125':'float32', 'AC_250':'float32', 'AC_500':'float32',
'AC_1000': 'float32', 'AC_2000': 'float32',
'AC_4000': 'float32', 'AC_8000': 'float32'}
df = df.astype(type_dict)
df = df.drop(columns='Study_ID')
df = df.fillna(pd.NA)
return df
def select_bp110(df):
# select BP110 data
is_bp110 = df['Device']=='BP110'
df_bp110 = df[is_bp110]
df_bp110.pop('Device')
return(df_bp110)
def select_bh5(df):
# select BAHA5P data
is_baha5p = df['Device']=='BAHA5P'
df_baha5p = df[is_baha5p]
df_baha5p.pop('Device')
return(df_baha5p)
act = read_AC_thresholds()
dsl = read_diff_sl()
act1 = act[['Device', 'AC_1000', 'AC_2000', 'AC_4000']].copy()
dsl1 = dsl[['1000_Hz', '2000_Hz', '4000_Hz']].copy()
combi = pd.concat([act1, dsl1], axis=1)
combi = combi.dropna()
combi.head()
Device | AC_1000 | AC_2000 | AC_4000 | 1000_Hz | 2000_Hz | 4000_Hz | |
---|---|---|---|---|---|---|---|
0 | BP110 | 10.0 | 10.0 | 0.0 | -15 | -15 | -12.0 |
1 | BP110 | 10.0 | 15.0 | 65.0 | -9 | 10 | 14.0 |
2 | BP110 | 15.0 | 10.0 | 30.0 | 7 | 6 | -16.0 |
3 | BP110 | 5.0 | 20.0 | 45.0 | -6 | 12 | 21.0 |
4 | BP110 | 5.0 | 5.0 | 15.0 | 6 | 6 | 13.0 |
combi_bh5 = select_bh5(combi)
combi_bp110 = select_bp110(combi)
combi_bh5.head()
AC_1000 | AC_2000 | AC_4000 | 1000_Hz | 2000_Hz | 4000_Hz | |
---|---|---|---|---|---|---|
20 | 10.0 | 15.0 | 25.0 | -10 | 10 | 19.0 |
37 | 20.0 | 5.0 | 15.0 | -14 | 17 | -1.0 |
47 | 5.0 | 5.0 | 0.0 | -10 | 4 | 3.0 |
51 | 20.0 | 5.0 | 20.0 | 3 | 4 | 13.0 |
53 | 10.0 | 15.0 | 35.0 | -6 | 5 | 15.0 |
for i in range(1, 4):
ac = (combi.iloc[:,i]).to_numpy(dtype='float')
dsl = (combi.iloc[:,i+3]).to_numpy(dtype='float')
lrgrs = linregress(ac, dsl)
print(2**(i-1),'kHz', lrgrs)
1 kHz LinregressResult(slope=0.557858211187556, intercept=-12.47454382291355, rvalue=0.4425469864353523, pvalue=2.4972840824021178e-05, stderr=0.12483220265884717, intercept_stderr=1.652778452997554) 2 kHz LinregressResult(slope=0.46707366238553244, intercept=-0.8541726660139375, rvalue=0.3949566704044585, pvalue=0.00020084156656662207, stderr=0.1199782704827768, intercept_stderr=1.7403743156956406) 4 kHz LinregressResult(slope=0.4116677037958929, intercept=-0.7643123833229613, rvalue=0.5667506432644367, pvalue=1.910915443092634e-08, stderr=0.06608710888949483, intercept_stderr=2.020281740040052)
for i in range(3):
ac_bp110 = (combi_bp110.iloc[:,i]).to_numpy(dtype='float')
dsl_bp110 = (combi_bp110.iloc[:,i+3]).to_numpy(dtype='float')
lrgrs = linregress(ac_bp110, dsl_bp110)
print(2**i,'kHz', lrgrs)
1 kHz LinregressResult(slope=0.4683760683760683, intercept=-8.269230769230768, rvalue=0.38883911350472916, pvalue=0.006307663108920862, stderr=0.16362497082158356, intercept_stderr=2.240527187252648) 2 kHz LinregressResult(slope=0.5121700149427195, intercept=-4.33554706956666, rvalue=0.4360064303094715, pvalue=0.001950075679921986, stderr=0.15586823132274572, intercept_stderr=2.388885112674612) 4 kHz LinregressResult(slope=0.44133635334088306, intercept=-5.06703001132502, rvalue=0.6251011783809193, pvalue=2.0387955257368852e-06, stderr=0.0812526827341766, intercept_stderr=2.70630597533495)
for i in range(3):
ac_bh5 = (combi_bh5.iloc[:,i]).to_numpy(dtype='float')
dsl_bh5 = (combi_bh5.iloc[:,i+3]).to_numpy(dtype='float')
lrgrs = linregress(ac_bh5, dsl_bh5)
print(2**i,'kHz', lrgrs)
1 kHz LinregressResult(slope=0.6028816466552315, intercept=-17.195883361921098, rvalue=0.563547478083405, pvalue=0.00034505507478122747, stderr=0.1515605365795818, intercept_stderr=1.9112725801403203) 2 kHz LinregressResult(slope=0.5303977272727272, intercept=2.301136363636365, rvalue=0.5010167375266132, pvalue=0.0018553647145924182, stderr=0.1571251545670172, intercept_stderr=2.0950020608935627) 4 kHz LinregressResult(slope=0.48914498141263935, intercept=2.1301115241635706, rvalue=0.6571466845317796, pvalue=1.3385378867428729e-05, stderr=0.09622113497938183, intercept_stderr=2.5482962561189266)
ac_1k_bh5 = combi_bh5.iloc[:,0].to_numpy(dtype='float')
dsl_1k_bh5 = (combi_bh5.iloc[:,3]).to_numpy(dtype='float')
lrgrs_1k_bh5 = linregress(ac_1k_bh5, dsl_1k_bh5)
r2_1k_bh5 = (lrgrs_1k_bh5[2])**2
r2_1k_bh5 = round(r2_1k_bh5, 2)
s_1k = 'r square 1 kHz =' + str(r2_1k_bh5)
print(s_1k)
r square 1 kHz =0.32
plt.scatter(ac_1k_bh5, dsl_1k_bh5, c='black')
plt.title("Relation between AC thresholds and diff SL at 1 kHz")
plt.xlabel("AC thresholds at 1 kHz [dB HL]")
plt.ylabel("Diff. SL BC-AC at 1 kHz [dB]")
slope = lrgrs_1k_bh5[0]
intercept = lrgrs_1k_bh5[1]
plt.plot(ac_1k_bh5, slope*ac_1k_bh5 + intercept, color='grey')
plt.text(0, 0, s_1k)
plt.show()
ac_2k_bh5 = combi_bh5.iloc[:,1].to_numpy(dtype='float')
dsl_2k_bh5 = (combi_bh5.iloc[:,4]).to_numpy(dtype='float')
lrgrs_2k_bh5 = linregress(ac_2k_bh5, dsl_2k_bh5)
r2_2k_bh5= (lrgrs_2k_bh5[2])**2
r2_2k_bh5 = round(r2_2k_bh5, 2)
s_2k = 'r square 2 kHz =' + str(r2_2k_bh5)
print(s_2k)
r square 2 kHz =0.25
plt.scatter(ac_2k_bh5, dsl_2k_bh5, c='black')
plt.title("Relation between AC thresholds and diff SL at 2 kHz")
plt.xlabel("AC thresholds at 2 kHz [dB HL]")
plt.ylabel("Diff. SL BC-AC at 2 kHz [dB]")
slope = lrgrs_2k_bh5[0]
intercept = lrgrs_2k_bh5[1]
plt.plot(ac_2k_bh5, slope*ac_2k_bh5 + intercept, color='grey')
plt.text(17, 7, s_2k)
plt.show()
ac_4k_bh5 = combi_bh5.iloc[:,2].to_numpy(dtype='float')
dsl_4k_bh5 = (combi_bh5.iloc[:,5]).to_numpy(dtype='float')
lrgrs_4k_bh5 = linregress(ac_4k_bh5, dsl_4k_bh5)
r2_4k_bh5= (lrgrs_4k_bh5[2])**2
r2_4k_bh5 = round(r2_4k_bh5, 2)
s_4k = 'r square 4 kHz =' + str(r2_4k_bh5)
print(s_4k)
r square 4 kHz =0.43
plt.scatter(ac_4k_bh5, dsl_4k_bh5, c='black')
plt.title("Relation between AC thresholds and diff SL at 4 kHz")
plt.xlabel("AC thresholds at 4 kHz [dB HL]")
plt.ylabel("Diff. SL BC-AC at 4 kHz [dB]")
slope = lrgrs_4k_bh5[0]
intercept = lrgrs_4k_bh5[1]
plt.plot(ac_4k_bh5, slope*ac_4k_bh5 + intercept, color='grey')
plt.text(34, 12, s_4k)
plt.show()