From d26d12e7081001bb907353237f9d0a13381433c6 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Sat, 10 Aug 2024 22:10:52 +0000 Subject: [PATCH] #21 add itemised trace to hERG QC --- pcpostprocess/hergQC.py | 164 +++++++++++++++++++++++++++------------- tests/test_herg_qc.py | 7 +- 2 files changed, 117 insertions(+), 54 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index 049e010..719f206 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -123,38 +123,57 @@ def run_qc(self, voltage_steps, times, before = self.filter_capacitive_spikes(before, times, voltage_steps) after = self.filter_capacitive_spikes(after, times, voltage_steps) + QC = {label: [(False, None)] for label in self.qc_labels} + if len(before) == 0 or len(after) == 0: - return False, [False for lab in self.qc_labels] + return QC if (None in qc_vals_before) or (None in qc_vals_after): - return False, False * self.no_QC + return QC qc1_1 = self.qc1(*qc_vals_before) qc1_2 = self.qc1(*qc_vals_after) - qc1 = [i and j for i, j in zip(qc1_1, qc1_2)] - qc2_1 = True - qc2_2 = True + QC['qc1.rseal'] = [qc1_1[0], qc1_2[0]] + QC['qc1.cm'] = [qc1_1[1], qc1_2[1]] + QC['qc1.rseries'] = [qc1_1[2], qc1_2[2]] + + QC['qc2.raw'] = [] + QC['qc2.subtracted'] = [] for i in range(n_sweeps): - qc2_1 = qc2_1 and self.qc2(before[i]) - qc2_2 = qc2_2 and self.qc2(before[i] - after[i]) + qc2_1 = self.qc2(before[i]) + qc2_2 = self.qc2(before[i] - after[i]) + + QC['qc2.raw'].append(qc2_1) + QC['qc2.subtracted'].append(qc2_2) qc3_1 = self.qc3(before[0, :], before[1, :]) qc3_2 = self.qc3(after[0, :], after[1, :]) qc3_3 = self.qc3(before[0, :] - after[0, :], before[1, :] - after[1, :]) + QC['qc3.raw'] = [qc3_1] + QC['qc3.E4031'] = [qc3_2] + QC['qc3.subtracted'] = [qc3_3] + rseals = [qc_vals_before[0], qc_vals_after[0]] cms = [qc_vals_before[1], qc_vals_after[1]] rseriess = [qc_vals_before[2], qc_vals_after[2]] qc4 = self.qc4(rseals, cms, rseriess) + QC['qc4.rseal'] = [qc4[0]] + QC['qc4.cm'] = [qc4[1]] + QC['qc4.rseries'] = [qc4[2]] + # indices where hERG peaks qc5 = self.qc5(before[0, :], after[0, :], self.qc5_win) qc5_1 = self.qc5_1(before[0, :], after[0, :], label='1') + QC['qc5.staircase'] = [qc5] + QC['qc5.1.staircase'] = [qc5_1] + # Ensure thats the windows are correct by checking the voltage trace assert np.all( np.abs(self.voltage[self.qc6_win[0]: self.qc6_win[1]] - 40.0))\ @@ -166,14 +185,17 @@ def run_qc(self, voltage_steps, times, np.abs(self.voltage[self.qc6_2_win[0]: self.qc6_2_win[1]] - 40.0))\ < 1e-8 - qc6, qc6_1, qc6_2 = True, True, True + QC['qc6.subtracted'] = [] + QC['qc6.1.subtracted'] = [] + QC['qc6.2.subtracted'] = [] for i in range(before.shape[0]): - qc6 = qc6 and self.qc6((before[i, :] - after[i, :]), - self.qc6_win, label='0') - qc6_1 = qc6_1 and self.qc6((before[i, :] - after[i, :]), - self.qc6_1_win, label='1') - qc6_2 = qc6_2 and self.qc6((before[i, :] - after[i, :]), - self.qc6_2_win, label='2') + qc6 = self.qc6((before[i, :] - after[i, :]), self.qc6_win, label="0") + qc6_1 = self.qc6((before[i, :] - after[i, :]), self.qc6_1_win, label="1") + qc6_2 = self.qc6((before[i, :] - after[i, :]), self.qc6_2_win, label="2") + + QC['qc6.subtracted'].append(qc6) + QC['qc6.1.subtracted'].append(qc6_1) + QC['qc6.2.subtracted'].append(qc6_2) if self._debug: fig = plt.figure(figsize=(8, 5)) @@ -199,37 +221,47 @@ def run_qc(self, voltage_steps, times, fig.savefig(os.path.join(self.plot_dir, 'qc_debug.png')) plt.close(fig) - # Make a flat list of all QC criteria (pass/fail bool) - QC = np.hstack([qc1, [qc2_1, qc2_2, qc3_1, qc3_2, qc3_3], - qc4, [qc5, qc5_1, qc6, qc6_1, qc6_2]]).flatten() + # Check if all QC criteria passed + passed = all([x for qc in QC.values() for x, _ in qc]) - passed = np.all(QC) return passed, QC def qc1(self, rseal, cm, rseries): - - if any([x is None for x in (rseal, cm, rseries)]): - return False, False, False - # Check R_seal, C_m, R_series within desired range - if rseal < self.rsealc[0] or rseal > self.rsealc[1] \ - or not np.isfinite(rseal): + if ( + rseal is None + or rseal < self.rsealc[0] + or rseal > self.rsealc[1] + or not np.isfinite(rseal) + ): self.logger.debug(f"rseal: {rseal}") qc11 = False else: qc11 = True - if cm < self.cmc[0] or cm > self.cmc[1] or not np.isfinite(cm): + + if ( + cm is None + or cm < self.cmc[0] + or cm > self.cmc[1] + or not np.isfinite(cm) + ): self.logger.debug(f"cm: {cm}") qc12 = False else: qc12 = True - if rseries < self.rseriesc[0] or rseries > self.rseriesc[1] \ - or not np.isfinite(rseries): + + if ( + rseries is None + or rseries < self.rseriesc[0] + or rseries > self.rseriesc[1] + or not np.isfinite(rseries) + ): self.logger.debug(f"rseries: {rseries}") qc13 = False else: qc13 = True - return [qc11, qc12, qc13] + + return [(qc11, rseal), (qc12, cm), (qc13, rseries)] def qc2(self, recording, method=3): # Check SNR is good @@ -245,9 +277,11 @@ def qc2(self, recording, method=3): if snr < self.snrc or not np.isfinite(snr): self.logger.debug(f"snr: {snr}") - return False + result = False + else: + result = True - return True + return (result, snr) def qc3(self, recording1, recording2, method=3): # Check 2 sweeps similar @@ -270,36 +304,51 @@ def qc3(self, recording1, recording2, method=3): rmsd = np.sqrt(np.mean((recording1 - recording2) ** 2)) if rmsd > rmsdc or not (np.isfinite(rmsd) and np.isfinite(rmsdc)): self.logger.debug(f"rmsd: {rmsd}, rmsdc: {rmsdc}") - return False - return True - - def qc4(self, rseals, cms, rseriess): + result = False + else: + result = True - if any([x is None for x in list(rseals) + list(cms) + list(rseriess)]): - return [False, False, False] + return (result, rmsd) + def qc4(self, rseals, cms, rseriess): # Check R_seal, C_m, R_series stability # Require std/mean < x% qc41 = True qc42 = True qc43 = True - if np.std(rseals) / np.mean(rseals) > self.rsealsc or not ( - np.isfinite(np.mean(rseals)) and np.isfinite(np.std(rseals))): - self.logger.debug(f"d_rseal: {np.std(rseals) / np.mean(rseals)}") + + if None in list(rseals): qc41 = False + d_rseal = None + else: + d_rseal = np.std(rseals) / np.mean(rseals) + if d_rseal > self.rsealsc or not ( + np.isfinite(np.mean(rseals)) and np.isfinite(np.std(rseals))): + self.logger.debug(f"d_rseal: {d_rseal}") + qc41 = False - if np.std(cms) / np.mean(cms) > self.cmsc or not ( - np.isfinite(np.mean(cms)) and np.isfinite(np.std(cms))): - self.logger.debug(f"d_cm: {np.std(cms) / np.mean(cms)}") + if None in list(cms): qc42 = False + d_cm = None + else: + d_cm = np.std(cms) / np.mean(cms) + if d_cm > self.cmsc or not ( + np.isfinite(np.mean(cms)) and np.isfinite(np.std(cms))): + self.logger.debug(f"d_cm: {d_cm}") + qc42 = False - if np.std(rseriess) / np.mean(rseriess) > self.rseriessc or not ( - np.isfinite(np.mean(rseriess)) - and np.isfinite(np.std(rseriess))): - self.logger.debug(f"d_rseries: {np.std(rseriess) / np.mean(rseriess)}") + if None in list(rseriess): qc43 = False + d_rseries = None + else: + d_rseries = np.std(rseriess) / np.mean(rseriess) + if d_rseries > self.rseriessc or not ( + np.isfinite(np.mean(rseriess)) + and np.isfinite(np.std(rseriess))): + self.logger.debug(f"d_rseries: {d_rseries}") + qc43 = False - return [qc41, qc42, qc43] + return [(qc41, d_rseal), (qc42, d_cm), (qc43, d_rseries)] def qc5(self, recording1, recording2, win=None, label=''): # Check pharma peak value drops after E-4031 application @@ -325,8 +374,11 @@ def qc5(self, recording1, recording2, win=None, label=''): if (max_diff < max_diffc) or not (np.isfinite(max_diff) and np.isfinite(max_diffc)): self.logger.debug(f"max_diff: {max_diff}, max_diffc: {max_diffc}") - return False - return True + result = False + else: + result = True + + return (result, max_diff) def qc5_1(self, recording1, recording2, win=None, label=''): # Check RMSD_0 drops after E-4031 application @@ -355,8 +407,11 @@ def qc5_1(self, recording1, recording2, win=None, label=''): if (rmsd0_diff < rmsd0_diffc) or not (np.isfinite(rmsd0_diff) and np.isfinite(rmsd0_diffc)): self.logger.debug(f"rmsd0_diff: {rmsd0_diff}, rmsd0c: {rmsd0_diffc}") - return False - return True + result = False + else: + result = True + + return (result, rmsd0_diff) def qc6(self, recording1, win=None, label=''): # Check subtracted staircase +40mV step up is non negative @@ -377,8 +432,11 @@ def qc6(self, recording1, win=None, label=''): if (val < valc) or not (np.isfinite(val) and np.isfinite(valc)): self.logger.debug(f"qc6_{label} val: {val}, valc: {valc}") - return False - return True + result = False + else: + result = True + + return (result, val) def filter_capacitive_spikes(self, current, times, voltage_step_times): """ diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 7d5bea9..f82b0f7 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -89,4 +89,9 @@ def test_run_qc(self): qc_vals_after_well, n_sweeps=2) logging.debug(well, passed) - self.assertTrue(passed) + + trace = "" + for label, results in qcs.items(): + if any([x == False for x, _ in results]): + trace += f"{label}: {results}\n" + self.assertTrue(passed, trace)