25 January 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
def read_AC_thresholds():
# read measured thresholds
f_in = '/media/guido/LACIE/Cingle_Guido/Master/AC_thresholds.xlsx'
p_in = Path(f_in)
df = pd.read_excel(p_in, sheet_name=0, header=0, index_col= 1, nrows=85)
df = df.drop(['Unnamed: 0'], axis=1)
type_dict = {'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)
new_dict ={'AC_500': 'AC_500_Hz', 'AC_1000': 'AC_1000_Hz',
'AC_2000': 'AC_2000_Hz', 'AC_4000': 'AC_4000_Hz'}
df = df.rename(columns=new_dict)
df = df.fillna(pd.NA)
return df
def present_AC_thresholds():
ac = read_AC_thresholds()
ac1 = ac[['Device', 'AC_500_Hz', 'AC_1000_Hz',
'AC_2000_Hz', 'AC_4000_Hz']].copy()
pta = ac1.iloc[:, 1:].mean(axis='columns')
ac1.insert(5, column= 'PTA (0.5-4 kHz)', value=pta)
pta_hf = ac1.iloc[:, 2:].mean(axis='columns')
ac1.insert(6, column= 'PTA_HF (1-4 kHz)', value=pta_hf)
return ac1
def read_BC_direct():
f_in = '/media/guido/LACIE/Cingle_Guido/Master/BC_direct_thresholds.xlsx'
p_in = Path(f_in)
df = pd.read_excel(p_in, sheet_name=0, header=0, index_col= 1, nrows=85)
df = df.drop(['Unnamed: 0'], axis=1)
type_dict = {'Device': 'str',
'BCd_250':'float32', 'BCd_250':'float32', 'BCd_500':'float32',
'BCd_1000': 'float32', 'BCd_1500': 'float32', 'BCd_2000': 'float32', 'BCd_3000': 'float32',
'BCd_4000': 'float32'}
new_dict ={'BCd_500': 'BCd_500_Hz', 'BCd_1000': 'BCd_1000_Hz',
'BCd_2000': 'BCd_2000_Hz', 'BCd_4000': 'BCd_4000_Hz'}
df = df.astype(type_dict)
df = df.rename(columns=new_dict)
df = df.fillna(pd.NA)
return df
def present_BCdirect_thresholds():
bc = read_BC_direct()
bc1 = bc[['Device', 'BCd_500_Hz', 'BCd_1000_Hz', 'BCd_2000_Hz', 'BCd_4000_Hz']].copy()
pta = bc1.iloc[:, 1:].mean(axis='columns')
bc1.insert(5, column= 'PTA (0.5-4 kHz)', value=pta)
pta_hf = bc1.iloc[:, 2:].mean(axis='columns')
bc1.insert(6, column= 'PTA_HF (1-4 kHz)', value=pta_hf)
bc1.dropna(inplace=True)
return bc1
def read_SII():
f_in = '/media/guido/LACIE/Cingle_Guido/Master/SII_Sfront65_Nfront65.xlsx'
p_in = Path(f_in)
df = pd.read_excel(p_in, header=0, index_col= 1, nrows=85)
df = df.drop(['Unnamed: 0',], axis=1)
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)
def select_cases():
sii = read_SII()
idx = sii.index
return idx
def select_ac_thresholds():
act = present_AC_thresholds()
idx = select_cases()
act1 = act.loc[idx]
return act1
def select_bc_thresholds():
bct = present_BCdirect_thresholds()
idx = select_cases()
bct1 = bct.loc[idx]
return bct1
# SII combi per Device type, as numpy array, presented as tuple
def sii_device():
sii = read_SII()
sii_bp110 = (select_bp110(sii).iloc[:, 2]).to_numpy(dtype='float')
sii_bh5 = (select_bh5(sii).iloc[:, 2]).to_numpy(dtype='float')
t = (sii_bp110, sii_bh5)
return t
# AC thresholds per Device type, as numpy array, presented as tuple
def act_device():
act = select_ac_thresholds()
act_bp110 = select_bp110(act)
act_bh5 = select_bh5(act)
t = (act_bp110, act_bh5)
return t
# BC direct thresholds per Device type, as numpy array, presented as tuple
def bct_device():
bct = select_bc_thresholds()
bct_bp110 = select_bp110(bct)
bct_bh5 = select_bh5(bct)
t = (bct_bp110, bct_bh5)
return t
# lin. reg. SII combination path (si) function thresholds (df)
def lreg_SII(thr, si):
l = list(thr.columns)
d = dict()
for j in range(6):
thres = (thr.iloc[:, j]).to_numpy(dtype='float')
lrgrs = linregress(thres, si)
r = lrgrs[2]
r2 = round(r**2, 3)
r = round(r, 3)
p = round(lrgrs[3], 3)
se = round(lrgrs[4], 4)
print(l[j], 'r square =', r2, ' r =', r, ' p =', p, ' se =', se)
de = {l[j] : (r, r2)}
d.update(de)
df = pd.DataFrame.from_dict(d, orient='index')
return df
(sii_bp110, sii_bh5) = sii_device()
(act_bp110, act_bh5) = act_device()
(bct_bp110, bct_bh5) = bct_device()
# BAHA5P: relation sii combination path vs. AC thresholds
lr_bh5_ac = lreg_SII(act_bh5, sii_bh5)
AC_500_Hz r square = 0.165 r = -0.407 p = 0.015 se = 0.0001 AC_1000_Hz r square = 0.084 r = -0.29 p = 0.091 se = 0.0001 AC_2000_Hz r square = 0.178 r = -0.422 p = 0.012 se = 0.0001 AC_4000_Hz r square = 0.404 r = -0.635 p = 0.0 se = 0.0 PTA (0.5-4 kHz) r square = 0.39 r = -0.625 p = 0.0 se = 0.0001 PTA_HF (1-4 kHz) r square = 0.391 r = -0.625 p = 0.0 se = 0.0001
# BP110: relation sii combination path vs. AC thresholds
lr_bp110_ac = lreg_SII(act_bp110, sii_bp110)
AC_500_Hz r square = 0.005 r = 0.069 p = 0.693 se = 0.0001 AC_1000_Hz r square = 0.012 r = -0.111 p = 0.527 se = 0.0001 AC_2000_Hz r square = 0.122 r = -0.35 p = 0.039 se = 0.0001 AC_4000_Hz r square = 0.515 r = -0.718 p = 0.0 se = 0.0 PTA (0.5-4 kHz) r square = 0.351 r = -0.593 p = 0.0 se = 0.0001 PTA_HF (1-4 kHz) r square = 0.449 r = -0.67 p = 0.0 se = 0.0001
# BAHA5P: relation sii combination path vs. BC direct thresholds
lr_bh5_bc = lreg_SII(bct_bh5, sii_bh5)
BCd_500_Hz r square = 0.043 r = -0.208 p = 0.23 se = 0.0001 BCd_1000_Hz r square = 0.008 r = 0.088 p = 0.616 se = 0.0001 BCd_2000_Hz r square = 0.255 r = -0.505 p = 0.002 se = 0.0001 BCd_4000_Hz r square = 0.364 r = -0.603 p = 0.0 se = 0.0 PTA (0.5-4 kHz) r square = 0.249 r = -0.499 p = 0.002 se = 0.0001 PTA_HF (1-4 kHz) r square = 0.285 r = -0.534 p = 0.001 se = 0.0001
# BP110: relation sii combination path vs. BC direct thresholds
lr_bp110_bc = lreg_SII(bct_bp110, sii_bp110)
BCd_500_Hz r square = 0.0 r = -0.019 p = 0.915 se = 0.0001 BCd_1000_Hz r square = 0.001 r = -0.027 p = 0.879 se = 0.0001 BCd_2000_Hz r square = 0.147 r = -0.383 p = 0.023 se = 0.0001 BCd_4000_Hz r square = 0.609 r = -0.781 p = 0.0 se = 0.0 PTA (0.5-4 kHz) r square = 0.323 r = -0.569 p = 0.0 se = 0.0001 PTA_HF (1-4 kHz) r square = 0.428 r = -0.654 p = 0.0 se = 0.0001
pta4_bh5 = (bct_bh5.iloc[:, 5]).to_numpy(dtype='float')
lrgrs = linregress(pta4_bh5, sii_bh5)
slope = lrgrs[0]
intercept = lrgrs[1]
r2 = round(lrgrs[2]**2, 3)
s = 'r square = '+ str(r2)
plt.scatter(pta4_bh5, sii_bh5, c='black')
plt.title("BAHA5P: Relation between BC direct thresholds PTA (0.5-4 kHz) and SII combination path")
plt.xlabel("BC direct thresholds PTA (1-4 kHz [dB]")
plt.ylabel("SII combination path")
plt.plot(pta4_bh5, slope*pta4_bh5 + intercept, color='grey')
plt.text(42, 0.482, s)
plt.show()
ac4_bh5 = (act_bh5.iloc[:, 5]).to_numpy(dtype='float')
lrgrs = linregress(ac4_bh5, sii_bh5)
slope = lrgrs[0]
intercept = lrgrs[1]
r2 = round(lrgrs[2]**2, 3)
s = 'r square = '+ str(r2)
plt.scatter(ac4_bh5, sii_bh5, c='black')
plt.title("BAHA5P: Relation between AC PTA (0.5-4 kHz) and SII combination path")
plt.xlabel("BC direct thresholds PTA (0.5-4 kHz [dB]")
plt.ylabel("SII combination path")
plt.plot(ac4_bh5, slope*ac4_bh5 + intercept, color='grey')
plt.text(19, 0.487, s)
plt.show()