From 57ac313403ec21845bbfa7563ee7baea72b95e2b Mon Sep 17 00:00:00 2001
From: Layla Riemann <layla.riemann@ptb.de>
Date: Mon, 2 Aug 2021 09:33:24 +0000
Subject: [PATCH] python code BA plots
---
bland_altman.py | 123 ++++++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 123 insertions(+)
create mode 100644 bland_altman.py
diff --git a/bland_altman.py b/bland_altman.py
new file mode 100644
index 0000000..a83679f
--- /dev/null
+++ b/bland_altman.py
@@ -0,0 +1,123 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Mon Mar 16 15:42:58 2020
+
+@author: rieman01
+
+Code to obtain the Bland-Altman plots for the spectral shape for the HS-, GOIA-,
+and WURST-SPECIAL sequence variants.
+"""
+
+import matplotlib.pyplot as plt
+import numpy as np
+import glob
+import matplotlib
+from xlwt import Workbook
+
+font = {'family' : 'normal',
+ 'weight' : 'normal',
+ 'size' : 15}
+
+matplotlib.rc('font', **font)
+plt.close('all')
+
+def bland_altman_plot(data1, data2, a,b,col,coli, *args, **kwargs):
+ """to generate points in BA plots for the differnet pulses and scenarios;
+ blue = HS, orange = GOIA, green = WURST, upper row: R_0, middle row: R_1,Mc,
+ bottom row: R_1,Wc """
+ data1 = np.asarray(data1)
+ data2 = np.asarray(data2)
+ meani = np.mean(np.sum([abs(data1), abs(data2)], axis=1))
+ diffabs=np.sum(abs(data1))-np.sum(abs(data2))
+ axs[a,b].scatter(meani, diffabs, color=col,marker='o',edgecolor=coli, *args, **kwargs)
+ return diffabs
+
+def alterman(diff,a,b):
+ """to generate mean and +/- 1.96*SD line in plots"""
+ md = np.mean(diff) # Mean of the difference
+ sd = np.std(diff, ddof=1) # Standard deviation of the difference
+ print(np.round(md,3),np.round(sd,4))
+ axs[a,b].axhline(md, color='red', linestyle='--')
+ axs[a,b].axhline(md + (1.96*sd), color='gray', linestyle='--')
+ axs[a,b].axhline(md - (1.96*sd), color='gray', linestyle='--')
+
+files= glob.glob('*.npy')
+data=[]
+for i in files:
+ data.append(np.load(i))
+data=np.asarray(data)
+data=np.reshape(data,[int(len(data)/12),3,4,np.shape(data)[1]])
+fig, axs = plt.subplots(3, 3,sharex='all', sharey='all')
+HS_same=[]
+HS_repro=[]
+HS_2sess=[]
+GOIA_same=[]
+GOIA_repro=[]
+GOIA_2sess=[]
+WURST_same=[]
+WURST_repro=[]
+WURST_2sess=[]
+RM_same=[]
+RM_repro=[]
+RM_2sess=[]
+
+colors=['blue','orange','green']
+colori=['deepskyblue','orangered','limegreen']
+for i in range(len(data)):
+ for j in range(3):
+ diff_same=bland_altman_plot(data[i][j][2], data[i][j][3],0,j,colors[j],colori[j]) #R_0
+ diff_2sess3=bland_altman_plot(data[i][j][1], data[i][j][3],2,j,colors[j],colori[j]) #R_1,Mc
+ diff_repro=bland_altman_plot(data[i][j][0], data[i][j][1],1,j,colors[j],colori[j]) #R_1,Wc
+ if j==1:
+ GOIA_same.append(diff_same)
+ GOIA_repro.append(diff_repro)
+ GOIA_2sess.append(diff_2sess3)
+ elif j==0:
+ HS_same.append(diff_same)
+ HS_repro.append(diff_repro)
+ HS_2sess.append(diff_2sess3)
+ elif j==2:
+ WURST_same.append(diff_same)
+ WURST_repro.append(diff_repro)
+ WURST_2sess.append(diff_2sess3)
+
+for k in range(3):
+ if k==1:
+ alterman(GOIA_same,0,k)
+ alterman(GOIA_repro,1,k)
+ alterman(GOIA_2sess,2,k)
+ elif k==0:
+ alterman(HS_same,0,k)
+ alterman(HS_repro,1,k)
+ alterman(HS_2sess,2,k)
+ elif k==2:
+ alterman(WURST_same,0,k)
+ alterman(WURST_repro,1,k)
+ alterman(WURST_2sess,2,k)
+
+fig.text(0.03, 0.55, '$BA_{i,y}$ / a.u.', ha='center', va='center', rotation='vertical')
+fig.text(0.55, 0.019, '$BA_{i,x}$ / a.u.', ha='center', va='center')
+plt.subplots_adjust(top=0.987,bottom=0.14,left=0.14,right=0.987, wspace=0.05, hspace=0.05)
+plt.show()
+
+data_neu=np.sum(np.absolute(data),axis=3) #to generate .xls file for the REML analysis
+dat=['HS','GOIA','WURST']
+rep=['ASR','ESR','WBSR']
+sess=['1_1','1_2','2_1','2_2']
+wb = Workbook()
+sheet1 = wb.add_sheet('Sheet 1')
+sheet1.write(0,0,'Pulse')
+sheet1.write(0,1,'Subject')
+sheet1.write(0,2, 'Session')
+sheet1.write(0,3,'Data')
+nvol=9
+kk=0
+for i in range(3):
+ sheet1.write(kk+1,0, dat[i])
+ for j in range(nvol):
+ sheet1.write(kk+1,1, j+1)
+ for k in range(4):
+ sheet1.write(kk+1,2, sess[k])
+ sheet1.write(kk+1,3, data_neu[j,i,k])
+ kk+=1
+wb.save('pulse_shape.xls')
--
GitLab