From 264c69752418e15e9c451c7dfe15983e2dd01606 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 8 Aug 2024 14:02:14 +0100 Subject: [PATCH 01/33] #21 Fix qc_vals_after in herg qc test --- tests/test_herg_qc.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 9d0fb29..02ef35c 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -59,7 +59,7 @@ def test_run_qc(self): before = tr_before.get_trace_sweeps(sweeps) after = tr_after.get_trace_sweeps(sweeps) qc_vals_before = tr_before.get_onboard_QC_values(sweeps=sweeps) - qc_vals_after = tr_before.get_onboard_QC_values(sweeps=sweeps) + qc_vals_after = tr_after.get_onboard_QC_values(sweeps=sweeps) res = {} From b9ebc168db7f0fc7e894ed5d15a32d8baea9bd63 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 8 Aug 2024 16:09:49 +0000 Subject: [PATCH 02/33] #21 parameterise herg QC test --- tests/test_herg_qc.py | 93 ++++++++++++++++++++----------------------- 1 file changed, 43 insertions(+), 50 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 02ef35c..ede36d8 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -8,8 +8,11 @@ from pcpostprocess.hergQC import hERGQC +WELLS = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] + class TestHergQC(unittest.TestCase): + def setUp(self): filepath = os.path.join('tests', 'test_data', '13112023_MW2_FF', 'staircaseramp (2)_2kHz_15.01.07') @@ -32,68 +35,58 @@ def setUp(self): self.test_trace_after = Trace(filepath2, json_file2) def test_run_qc(self): - tr_before = self.test_trace_before - tr_after = self.test_trace_after - - v = tr_before.get_voltage() - - times = tr_after.get_times() - - self.assertTrue(np.all(np.isfinite(v))) - self.assertTrue(np.all(np.isfinite(times))) + for well in WELLS: + with self.subTest(well): + tr_before = self.test_trace_before + tr_after = self.test_trace_after - # Calculate sampling rate in (use kHz) - sampling_rate = int(1.0 / (times[1] - times[0])) + v = tr_before.get_voltage() - plot_dir = os.path.join(self.output_dir, - 'test_run_qc') + times = tr_after.get_times() - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) + self.assertTrue(np.all(np.isfinite(v))) + self.assertTrue(np.all(np.isfinite(times))) - hergqc = hERGQC(sampling_rate=sampling_rate, - plot_dir=plot_dir, - voltage=v) + # Calculate sampling rate in (use kHz) + sampling_rate = int(1.0 / (times[1] - times[0])) - sweeps = [0, 1] - before = tr_before.get_trace_sweeps(sweeps) - after = tr_after.get_trace_sweeps(sweeps) - qc_vals_before = tr_before.get_onboard_QC_values(sweeps=sweeps) - qc_vals_after = tr_after.get_onboard_QC_values(sweeps=sweeps) + plot_dir = os.path.join(self.output_dir, + 'test_run_qc') - res = {} + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) - # Spot check a few wells - # We could check all of the wells but it's time consuming + hergqc = hERGQC(sampling_rate=sampling_rate, + plot_dir=plot_dir, + voltage=v) - test_wells = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] + sweeps = [0, 1] + before = tr_before.get_trace_sweeps(sweeps) + after = tr_after.get_trace_sweeps(sweeps) + qc_vals_before = tr_before.get_onboard_QC_values(sweeps=sweeps) + qc_vals_after = tr_after.get_onboard_QC_values(sweeps=sweeps) - voltage_protocol = tr_before.get_voltage_protocol() + # Spot check a few wells + # We could check all of the wells but it's time consuming - for well in test_wells: - # Take values from the first sweep only - qc_vals_before_well = np.array(qc_vals_before[well])[0, :] - qc_vals_after_well = np.array(qc_vals_after[well])[0, :] + voltage_protocol = tr_before.get_voltage_protocol() - before_well = np.array(before[well]) - after_well = np.array(after[well]) + # Take values from the first sweep only + qc_vals_before_well = np.array(qc_vals_before[well])[0, :] + qc_vals_after_well = np.array(qc_vals_after[well])[0, :] - #  Assume that there are no discontinuities at the start or end of ramps - voltage_steps = [tstart - for tstart, tend, vstart, vend in - voltage_protocol.get_all_sections() if vend == vstart] + before_well = np.array(before[well]) + after_well = np.array(after[well]) - passed, qcs = hergqc.run_qc(voltage_steps, - times, before_well, after_well, - qc_vals_before_well, - qc_vals_after_well, n_sweeps=2) + #  Assume that there are no discontinuities at the start or end of ramps + voltage_steps = [tstart + for tstart, tend, vstart, vend in + voltage_protocol.get_all_sections() if vend == vstart] - logging.debug(well, passed) - res[well] = passed + passed, qcs = hergqc.run_qc(voltage_steps, + times, before_well, after_well, + qc_vals_before_well, + qc_vals_after_well, n_sweeps=2) - self.assertTrue(res['A01']) - self.assertTrue(res['A02']) - self.assertTrue(res['A03']) - self.assertTrue(res['A04']) - self.assertFalse(res['A05']) - self.assertFalse(res['D01']) + logging.debug(well, passed) + self.assertTrue(passed) From 3815c47e7fe78a1db231a20ba831397be4e0f8d1 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 8 Aug 2024 16:25:34 +0000 Subject: [PATCH 03/33] #21 simplify herg QC test parametrisation --- tests/test_herg_qc.py | 52 +++++++++++++++++++++---------------------- 1 file changed, 26 insertions(+), 26 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index ede36d8..7d5bea9 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -8,8 +8,6 @@ from pcpostprocess.hergQC import hERGQC -WELLS = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] - class TestHergQC(unittest.TestCase): @@ -35,40 +33,42 @@ def setUp(self): self.test_trace_after = Trace(filepath2, json_file2) def test_run_qc(self): - for well in WELLS: - with self.subTest(well): - tr_before = self.test_trace_before - tr_after = self.test_trace_after + tr_before = self.test_trace_before + tr_after = self.test_trace_after - v = tr_before.get_voltage() + v = tr_before.get_voltage() - times = tr_after.get_times() + times = tr_after.get_times() - self.assertTrue(np.all(np.isfinite(v))) - self.assertTrue(np.all(np.isfinite(times))) + self.assertTrue(np.all(np.isfinite(v))) + self.assertTrue(np.all(np.isfinite(times))) - # Calculate sampling rate in (use kHz) - sampling_rate = int(1.0 / (times[1] - times[0])) + # Calculate sampling rate in (use kHz) + sampling_rate = int(1.0 / (times[1] - times[0])) - plot_dir = os.path.join(self.output_dir, - 'test_run_qc') + plot_dir = os.path.join(self.output_dir, + 'test_run_qc') - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) - hergqc = hERGQC(sampling_rate=sampling_rate, - plot_dir=plot_dir, - voltage=v) + hergqc = hERGQC(sampling_rate=sampling_rate, + plot_dir=plot_dir, + voltage=v) - sweeps = [0, 1] - before = tr_before.get_trace_sweeps(sweeps) - after = tr_after.get_trace_sweeps(sweeps) - qc_vals_before = tr_before.get_onboard_QC_values(sweeps=sweeps) - qc_vals_after = tr_after.get_onboard_QC_values(sweeps=sweeps) + sweeps = [0, 1] + before = tr_before.get_trace_sweeps(sweeps) + after = tr_after.get_trace_sweeps(sweeps) + qc_vals_before = tr_before.get_onboard_QC_values(sweeps=sweeps) + qc_vals_after = tr_after.get_onboard_QC_values(sweeps=sweeps) - # Spot check a few wells - # We could check all of the wells but it's time consuming + # Spot check a few wells + # We could check all of the wells but it's time consuming + test_wells = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] + + for well in test_wells: + with self.subTest(well): voltage_protocol = tr_before.get_voltage_protocol() # Take values from the first sweep only From d26d12e7081001bb907353237f9d0a13381433c6 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Sat, 10 Aug 2024 22:10:52 +0000 Subject: [PATCH 04/33] #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) From 2f564791b54a3d4da2f59f646f61c45b23a1c04f Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 12 Aug 2024 10:59:34 +0000 Subject: [PATCH 05/33] #21 Annotate converted QC units --- pcpostprocess/hergQC.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index 719f206..a4cd154 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -38,9 +38,9 @@ def __init__(self, sampling_rate=5, plot_dir=None, voltage=np.array([]), # Define all thresholds # qc1 - self.rsealc = [1e8, 1e12] # in Ohm # TODO double check values - self.cmc = [1e-12, 1e-10] # in F - self.rseriesc = [1e6, 2.5e7] # in Ohm + self.rsealc = [1e8, 1e12] # in Ohm, converted from [0.1, 1000] GOhm + self.cmc = [1e-12, 1e-10] # in F, converted from [1, 100] pF + self.rseriesc = [1e6, 2.5e7] # in Ohm, converted from [1, 25] MOhm # qc2 self.snrc = 25 # qc3 From fa8644d88142d50976539a933624ad3e25e3ec3e Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 12 Aug 2024 13:43:09 +0000 Subject: [PATCH 06/33] #21 Update calls to hergQC methods --- pcpostprocess/hergQC.py | 6 ++++-- pcpostprocess/scripts/run_herg_qc.py | 20 ++++++++++++-------- tests/test_herg_qc.py | 4 ++-- 3 files changed, 18 insertions(+), 12 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index a4cd154..04e52f5 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -1,6 +1,8 @@ import logging import os +from collections import OrderedDict + import matplotlib.pyplot as plt import numpy as np import scipy.stats @@ -123,7 +125,7 @@ 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} + QC = OrderedDict([(label, [(False, None)]) for label in self.qc_labels]) if len(before) == 0 or len(after) == 0: return QC @@ -213,7 +215,7 @@ def run_qc(self, voltage_steps, times, plt.ylabel('Current [pA]') # https://stackoverflow.com/a/13589144 - from collections import OrderedDict # fix legend labels + # fix legend labels handles, labels = plt.gca().get_legend_handles_labels() by_label = OrderedDict(zip(labels, handles)) fig.legend(by_label.values(), by_label.keys()) diff --git a/pcpostprocess/scripts/run_herg_qc.py b/pcpostprocess/scripts/run_herg_qc.py index 237c08e..d592452 100644 --- a/pcpostprocess/scripts/run_herg_qc.py +++ b/pcpostprocess/scripts/run_herg_qc.py @@ -701,18 +701,22 @@ def extract_protocol(readname, savename, time_strs, selected_wells, args): row_dict['QC6'] = hergqc.qc6(current, win=hergqc.qc6_win, - label='0') + label='0')[0] #  Assume there is only one sweep for all non-QC protocols rseal_before, cm_before, rseries_before = qc_before[well][0] rseal_after, cm_after, rseries_after = qc_after[well][0] - row_dict['QC1'] = all(list(hergqc.qc1(rseal_before, cm_before, rseries_before)) + - list(hergqc.qc1(rseal_after, cm_after, rseries_after))) + qc1_1 = hergqc.qc1(rseal_before, cm_before, rseries_before) + qc1_2 = hergqc.qc1(rseal_after, cm_after, rseries_after) - row_dict['QC4'] = all(hergqc.qc4([rseal_before, rseal_after], - [cm_before, cm_after], - [rseries_before, rseries_after])) + row_dict['QC1'] = all([x for x, _ in qc1_1 + qc1_2]) + + qc4 = hergqc.qc4([rseal_before, rseal_after], + [cm_before, cm_after], + [rseries_before, rseries_after]) + + row_dict['QC4'] = all([x for x, _ in qc4]) if args.output_traces: out_fname = os.path.join(traces_dir, @@ -935,7 +939,7 @@ def run_qc_for_protocol(readname, savename, time_strs, args): np.array(qc_before[well])[0, :], np.array(qc_after[well])[0, :], nsweeps) - df_rows.append([well] + list(QC)) + df_rows.append([well] + [all([x for x, _ in qc]) for qc in QC.values()]) if selected: selected_wells.append(well) @@ -1086,7 +1090,7 @@ def qc3_bookend(readname, savename, time_strs, args): last_processed[well], times, voltage_steps ).flatten() - passed = hergqc.qc3(trace1, trace2) + passed = hergqc.qc3(trace1, trace2)[0] res_dict[well] = passed diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index f82b0f7..6728135 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -67,10 +67,10 @@ def test_run_qc(self): test_wells = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] + voltage_protocol = tr_before.get_voltage_protocol() + for well in test_wells: with self.subTest(well): - voltage_protocol = tr_before.get_voltage_protocol() - # Take values from the first sweep only qc_vals_before_well = np.array(qc_vals_before[well])[0, :] qc_vals_after_well = np.array(qc_vals_after[well])[0, :] From 78e172e0e17b491c34c4cc9a9edb3e1678d9aead Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 12 Aug 2024 15:57:44 +0000 Subject: [PATCH 07/33] #21 Add unit test for herg qc1 --- tests/test_herg_qc.py | 97 +++++++++++++++++++++++++++++++++++-------- 1 file changed, 80 insertions(+), 17 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 6728135..b129ff3 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -29,22 +29,85 @@ def setUp(self): if not os.path.exists(self.output_dir): os.makedirs(self.output_dir) + self.test_trace_before = Trace(filepath, json_file) self.test_trace_after = Trace(filepath2, json_file2) - def test_run_qc(self): - tr_before = self.test_trace_before - tr_after = self.test_trace_after + self.voltage = self.test_trace_before.get_voltage() + self.times = self.test_trace_after.get_times() - v = tr_before.get_voltage() + # Calculate sampling rate in (use kHz) + self.sampling_rate = int(1.0 / (self.times[1] - self.times[0])) - times = tr_after.get_times() + def test_qc1(self): + def passed(result): + return all([x for x, _ in result]) - self.assertTrue(np.all(np.isfinite(v))) - self.assertTrue(np.all(np.isfinite(times))) + plot_dir = os.path.join(self.output_dir, "test_qc1") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) - # Calculate sampling rate in (use kHz) - sampling_rate = int(1.0 / (times[1] - times[0])) + hergqc = hERGQC(sampling_rate=self.sampling_rate, + plot_dir=plot_dir, + voltage=self.voltage) + + # qc1 checks that rseal, cm, rseries are within range + rseal_lo, rseal_hi = 1e8, 1e12 + rseal_mid = (rseal_lo + rseal_hi) / 2 + + cm_lo, cm_hi = 1e-12, 1e-10 + cm_mid = (cm_lo + cm_hi) / 2 + + rseries_lo, rseries_hi = 1e6, 2.5e7 + rseries_mid = (rseries_lo + rseries_hi) / 2 + + tol = 1e-3 + + test_matrix = [ + [(rseal_lo, cm_lo, rseries_lo), True], + [(rseal_mid, cm_mid, rseries_mid), True], + [(rseal_hi, cm_hi, rseries_hi), True], + [(rseal_lo - tol, cm_lo, rseries_lo), False], + [(rseal_lo, cm_lo - tol, rseries_lo), False], + [(rseal_lo, cm_lo, rseries_lo - tol), False], + [(rseal_hi + tol, cm_hi, rseries_hi), False], + [(rseal_hi, cm_hi + tol, rseries_hi), False], + [(rseal_hi, cm_hi, rseries_hi + tol), False], + [(np.inf, cm_mid, rseries_mid), False], + [(rseal_mid, np.inf, rseries_mid), False], + [(rseal_mid, cm_mid, np.inf), False], + [(None, cm_mid, rseries_mid), False], + [(rseal_mid, None, rseries_mid), False], + [(rseal_mid, cm_mid, None), False], + [(None, None, None), False], + [(0, 0, 0), False], + ] + + for (rseal, cm, rseries), expected in test_matrix: + self.assertEqual( + passed(hergqc.qc1(rseal, cm, rseries)), + expected, + f"({rseal}, {cm}, {rseries})", + ) + + def test_qc2(self): + pass + + def test_qc3(self): + pass + + def test_qc4(self): + pass + + def test_qc5(self): + pass + + def test_qc6(self): + pass + + def test_run_qc(self): + self.assertTrue(np.all(np.isfinite(self.voltage))) + self.assertTrue(np.all(np.isfinite(self.times))) plot_dir = os.path.join(self.output_dir, 'test_run_qc') @@ -52,22 +115,22 @@ def test_run_qc(self): if not os.path.exists(plot_dir): os.makedirs(plot_dir) - hergqc = hERGQC(sampling_rate=sampling_rate, + hergqc = hERGQC(sampling_rate=self.sampling_rate, plot_dir=plot_dir, - voltage=v) + voltage=self.voltage) sweeps = [0, 1] - before = tr_before.get_trace_sweeps(sweeps) - after = tr_after.get_trace_sweeps(sweeps) - qc_vals_before = tr_before.get_onboard_QC_values(sweeps=sweeps) - qc_vals_after = tr_after.get_onboard_QC_values(sweeps=sweeps) + before = self.test_trace_before.get_trace_sweeps(sweeps) + after = self.test_trace_after.get_trace_sweeps(sweeps) + qc_vals_before = self.test_trace_before.get_onboard_QC_values(sweeps=sweeps) + qc_vals_after = self.test_trace_after.get_onboard_QC_values(sweeps=sweeps) # Spot check a few wells # We could check all of the wells but it's time consuming test_wells = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] - voltage_protocol = tr_before.get_voltage_protocol() + voltage_protocol = self.test_trace_before.get_voltage_protocol() for well in test_wells: with self.subTest(well): @@ -84,7 +147,7 @@ def test_run_qc(self): voltage_protocol.get_all_sections() if vend == vstart] passed, qcs = hergqc.run_qc(voltage_steps, - times, before_well, after_well, + self.times, before_well, after_well, qc_vals_before_well, qc_vals_after_well, n_sweeps=2) From bfe9047a3467fc6d2d5b1e4fd61ec6269832286a Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 12 Aug 2024 16:39:45 +0000 Subject: [PATCH 08/33] #21 Adjust tolerances for herg qc1 unit test --- tests/test_herg_qc.py | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index b129ff3..a8c535b 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -61,18 +61,20 @@ def passed(result): rseries_lo, rseries_hi = 1e6, 2.5e7 rseries_mid = (rseries_lo + rseries_hi) / 2 - tol = 1e-3 + rseal_tol = 0.1 + cm_tol = 1e-13 + rseries_tol = 0.1 test_matrix = [ [(rseal_lo, cm_lo, rseries_lo), True], [(rseal_mid, cm_mid, rseries_mid), True], [(rseal_hi, cm_hi, rseries_hi), True], - [(rseal_lo - tol, cm_lo, rseries_lo), False], - [(rseal_lo, cm_lo - tol, rseries_lo), False], - [(rseal_lo, cm_lo, rseries_lo - tol), False], - [(rseal_hi + tol, cm_hi, rseries_hi), False], - [(rseal_hi, cm_hi + tol, rseries_hi), False], - [(rseal_hi, cm_hi, rseries_hi + tol), False], + [(rseal_lo - rseal_tol, cm_lo, rseries_lo), False], + [(rseal_lo, cm_lo - cm_tol, rseries_lo), False], + [(rseal_lo, cm_lo, rseries_lo - rseries_tol), False], + [(rseal_hi + rseal_tol, cm_hi, rseries_hi), False], + [(rseal_hi, cm_hi + cm_tol, rseries_hi), False], + [(rseal_hi, cm_hi, rseries_hi + rseries_tol), False], [(np.inf, cm_mid, rseries_mid), False], [(rseal_mid, np.inf, rseries_mid), False], [(rseal_mid, cm_mid, np.inf), False], From e5cfda3349dc2ab8ff82ccb452ba110c2d2991f8 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 12 Aug 2024 20:39:27 +0000 Subject: [PATCH 09/33] #21 Add QCDict --- pcpostprocess/hergQC.py | 84 ++++++++++++++++++++-------- pcpostprocess/scripts/run_herg_qc.py | 19 +++---- tests/test_herg_qc.py | 18 +++--- 3 files changed, 79 insertions(+), 42 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index deb12ba..418a5dc 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -8,16 +8,67 @@ import scipy.stats -class hERGQC(object): +class QCDict(object): + + labels = [ + "qc1.rseal", + "qc1.cm", + "qc1.rseries", + "qc2.raw", + "qc2.subtracted", + "qc3.raw", + "qc3.E4031", + "qc3.subtracted", + "qc4.rseal", + "qc4.cm", + "qc4.rseries", + "qc5.staircase", + "qc5.1.staircase", + "qc6.subtracted", + "qc6.1.subtracted", + "qc6.2.subtracted", + ] + + def __init__(self): + self._dict = OrderedDict([(label, [(False, None)]) for label in QCDict.labels]) + + def __str__(self): + return self._dict.__str__() + + def __repr__(self): + return self._dict.__repr__() + + def __getitem__(self, key): + return self._dict.__getitem__(key) + + def __setitem__(self, key, value): + if key not in QCDict.labels: + raise KeyError(f"Invalid QC key: {key}") + self._dict.__setitem__(key, value) + + def keys(self): + return self._dict.keys() + + def items(self): + return self._dict.items() + + def values(self): + return self._dict.values() + + def qc_passed(self, label): + """Return whether a single QC passed.""" + return all([x for x, _ in self._dict[label]]) + + def passed_list(self): + """Return a list of booleans indicating whether each QC passed.""" + return [self.qc_passed(label) for label in QCDict.labels] + + def all_passed(self): + """Return whether all QC passed.""" + return all(self.passed_list()) - QCnames = ['qc1.rseal', 'qc1.cm', 'qc1.rseries', - 'qc2.raw', 'qc2.subtracted', - 'qc3.raw', 'qc3.E4031', 'qc3.subtracted', - 'qc4.rseal', 'qc4.cm', 'qc4.rseries', - 'qc5.staircase', 'qc5.1.staircase', - 'qc6.subtracted', 'qc6.1.subtracted', 'qc6.2.subtracted'] - no_QC = len(QCnames) +class hERGQC(object): def __init__(self, sampling_rate=5, plot_dir=None, voltage=np.array([]), n_sweeps=None, removal_time=5): @@ -82,16 +133,6 @@ def __init__(self, sampling_rate=5, plot_dir=None, voltage=np.array([]), self._debug = True - self.qc_labels = ['qc1.rseal', 'qc1.cm', 'qc1.rseries', 'qc2.raw', - 'qc2.subtracted', 'qc3.raw', 'qc3.E4031', - 'qc3.subtracted', 'qc4.rseal', 'qc4.cm', - 'qc4.rseries', 'qc5.staircase', 'qc5.1.staircase', - 'qc6.subtracted', 'qc6.1.subtracted', - 'qc6.2.subtracted'] - - def get_qc_names(self): - return self.QCnames - def set_trace(self, before, after, qc_vals_before, qc_vals_after, n_sweeps): self._before = before @@ -125,7 +166,7 @@ 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 = OrderedDict([(label, [(False, None)]) for label in self.qc_labels]) + QC = QCDict() if len(before) == 0 or len(after) == 0: return QC @@ -223,10 +264,7 @@ def run_qc(self, voltage_steps, times, fig.savefig(os.path.join(self.plot_dir, 'qc_debug.png')) plt.close(fig) - # Check if all QC criteria passed - passed = all([x for qc in QC.values() for x, _ in qc]) - - return passed, QC + return QC def qc1(self, rseal, cm, rseries): # Check R_seal, C_m, R_series within desired range diff --git a/pcpostprocess/scripts/run_herg_qc.py b/pcpostprocess/scripts/run_herg_qc.py index c496c33..85cedf1 100644 --- a/pcpostprocess/scripts/run_herg_qc.py +++ b/pcpostprocess/scripts/run_herg_qc.py @@ -922,16 +922,15 @@ def run_qc_for_protocol(readname, savename, time_strs, args): voltage_protocol.get_all_sections() if vend == vstart] # Run QC with raw currents - _, QC = hergqc.run_qc(voltage_steps, times, - before_currents_corrected, - after_currents_corrected, - np.array(qc_before[well])[0, :], - np.array(qc_after[well])[0, :], nsweeps) - - passed_qc = [all([x for x, _ in qc]) for qc in QC.values()] - df_rows.append([well] + passed_qc) - - selected = np.all(QC) and not no_cell + QC = hergqc.run_qc(voltage_steps, times, + before_currents_corrected, + after_currents_corrected, + np.array(qc_before[well])[0, :], + np.array(qc_after[well])[0, :], nsweeps) + + df_rows.append([well] + QC.passed_list()) + + selected = QC.all_passed() and not no_cell if selected: selected_wells.append(well) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index a8c535b..f4b9a0a 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -148,15 +148,15 @@ def test_run_qc(self): for tstart, tend, vstart, vend in voltage_protocol.get_all_sections() if vend == vstart] - passed, qcs = hergqc.run_qc(voltage_steps, - self.times, before_well, after_well, - qc_vals_before_well, - qc_vals_after_well, n_sweeps=2) + QC = hergqc.run_qc(voltage_steps, + self.times, before_well, after_well, + qc_vals_before_well, + qc_vals_after_well, n_sweeps=2) - logging.debug(well, passed) + logging.debug(well, QC.all_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) + for label, result in QC.items(): + if not QC.qc_passed(label): + trace += f"{label}: {result}\n" + self.assertTrue(QC.all_passed(), trace) From 27a6a76b67ef1f1a573df2e5bbda03a576ae18e4 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 18 Nov 2024 14:43:04 +0000 Subject: [PATCH 10/33] #21 Add test_qc2 --- tests/test_herg_qc.py | 25 ++++++++++++++++++++++++- 1 file changed, 24 insertions(+), 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index f4b9a0a..6e09d50 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -93,7 +93,30 @@ def passed(result): ) def test_qc2(self): - pass + plot_dir = os.path.join(self.output_dir, "test_qc2") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC(sampling_rate=self.sampling_rate, + plot_dir=plot_dir, + voltage=self.voltage) + + # qc2 checks that raw and subtracted SNR are above a minimum threshold + recording = np.asarray([0, 1] * 500 + [0, 10] * 500) # snr = 70.75 + result = hergqc.qc2(recording) + self.assertTrue(result[0], f"({result[1]})") + + recording = np.asarray([0, 1] * 500 + [0, 6.03125] * 500) # snr = 25.02 + result = hergqc.qc2(recording) + self.assertTrue(result[0], f"({result[1]})") + + recording = np.asarray([0, 1] * 500 + [0, 6.015625] * 500) # snr = 24.88 + result = hergqc.qc2(recording) + self.assertFalse(result[0], f"({result[1]})") + + recording = np.asarray([0, 1] * 1000) # snr = 1.0 + result = hergqc.qc2(recording) + self.assertFalse(result[0], f"({result[1]})") def test_qc3(self): pass From bed87c13fdd92fdc88a0dd6b7547086e0befc462 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 18 Nov 2024 15:10:22 +0000 Subject: [PATCH 11/33] #21 Add test_qc3 --- tests/test_herg_qc.py | 25 ++++++++++++++++++++++++- 1 file changed, 24 insertions(+), 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 6e09d50..a5f2ae7 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -119,7 +119,30 @@ def test_qc2(self): self.assertFalse(result[0], f"({result[1]})") def test_qc3(self): - pass + plot_dir = os.path.join(self.output_dir, "test_qc3") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC(sampling_rate=self.sampling_rate, + plot_dir=plot_dir, + voltage=self.voltage) + + # qc3 checks that rmsd of two sweeps are similar + test_matrix = [ + (0, True), + (1, True), + (2, True), + (3, True), + (4, False), + (5, False), + (6, False), + ] + + for i, expected in test_matrix: + recording0 = np.asarray([0, 1] * 1000) + recording1 = np.asarray(recording0 + i) + result = hergqc.qc3(recording0, recording1) + self.assertEqual(result[0], expected, f"({result[1]})") def test_qc4(self): pass From 14fd05224466fc18dd3c0082df4b0ce5dd985043 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 18 Nov 2024 15:40:21 +0000 Subject: [PATCH 12/33] #21 Add test_qc4 --- tests/test_herg_qc.py | 37 ++++++++++++++++++++++++++++++++++++- 1 file changed, 36 insertions(+), 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index a5f2ae7..f317bc1 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -92,6 +92,8 @@ def passed(result): f"({rseal}, {cm}, {rseries})", ) + # TODO: Test on select data + def test_qc2(self): plot_dir = os.path.join(self.output_dir, "test_qc2") if not os.path.exists(plot_dir): @@ -118,6 +120,8 @@ def test_qc2(self): result = hergqc.qc2(recording) self.assertFalse(result[0], f"({result[1]})") + # TODO: Test on select data + def test_qc3(self): plot_dir = os.path.join(self.output_dir, "test_qc3") if not os.path.exists(plot_dir): @@ -144,13 +148,44 @@ def test_qc3(self): result = hergqc.qc3(recording0, recording1) self.assertEqual(result[0], expected, f"({result[1]})") + # TODO: Test on select data + def test_qc4(self): - pass + def passed(result): + return all([x for x, _ in result]) + + plot_dir = os.path.join(self.output_dir, "test_qc1") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC( + sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage + ) + + # qc4 checks that rseal, cm, rseries are similar before/after E-4031 change + test_matrix = [ + [([1, 1], [1, 1], [1, 1]), True], + [([1, 2], [1, 2], [1, 2]), True], + [([1, 3], [1, 3], [1, 3]), True], + [([1, 4], [1, 4], [1, 4]), False], + [([1, 5], [1, 5], [1, 5]), False], + ] + + for (rseals, cms, rseriess), expected in test_matrix: + self.assertEqual( + passed(hergqc.qc4(rseals, cms, rseriess)), + expected, + f"({rseals}, {cms}, {rseriess})", + ) + + # TODO: Test on select data def test_qc5(self): + # TODO: Test on select data pass def test_qc6(self): + # TODO: Test on select data pass def test_run_qc(self): From aac368326811a0d021544496e79a37f19f7decb9 Mon Sep 17 00:00:00 2001 From: Joseph Date: Mon, 18 Nov 2024 16:07:53 +0000 Subject: [PATCH 13/33] update test --- tests/test_herg_qc.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index f4b9a0a..a15bb46 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -130,7 +130,7 @@ def test_run_qc(self): # Spot check a few wells # We could check all of the wells but it's time consuming - test_wells = ['A01', 'A02', 'A03', 'A04', 'A05', 'D01'] + test_wells = ['A01', 'A02', 'A03'] voltage_protocol = self.test_trace_before.get_voltage_protocol() @@ -159,4 +159,5 @@ def test_run_qc(self): for label, result in QC.items(): if not QC.qc_passed(label): trace += f"{label}: {result}\n" + print(f"Testing Well, {well}") self.assertTrue(QC.all_passed(), trace) From e138f4694bc8e62eacefce27c7038e7ead9a3bd8 Mon Sep 17 00:00:00 2001 From: Joseph Date: Mon, 18 Nov 2024 16:24:03 +0000 Subject: [PATCH 14/33] Add simplified script --- .../scripts/export_staircase_data.py | 730 ++++++++++++++++++ 1 file changed, 730 insertions(+) create mode 100644 pcpostprocess/scripts/export_staircase_data.py diff --git a/pcpostprocess/scripts/export_staircase_data.py b/pcpostprocess/scripts/export_staircase_data.py new file mode 100644 index 0000000..ddf7b20 --- /dev/null +++ b/pcpostprocess/scripts/export_staircase_data.py @@ -0,0 +1,730 @@ +import argparse +import datetime +import importlib.util +import json +import logging +import multiprocessing +import os +import string +import subprocess +import sys + +import cycler +import matplotlib +import matplotlib.pyplot as plt +import numpy as np +import pandas as pd +import regex as re +import scipy +from syncropatch_export.trace import Trace as SyncropatchTrace +from syncropatch_export.voltage_protocols import VoltageProtocol + +from pcpostprocess.detect_ramp_bounds import detect_ramp_bounds +from pcpostprocess.hergQC import hERGQC +from pcpostprocess.infer_reversal import infer_reversal_potential +from pcpostprocess.leak_correct import fit_linear_leak, get_leak_corrected +from pcpostprocess.subtraction_plots import do_subtraction_plot + +pool_kws = {'maxtasksperchild': 1} + +color_cycle = ["#5790fc", "#f89c20", "#e42536", "#964a8b", "#9c9ca1", "#7a21dd"] +plt.rcParams['axes.prop_cycle'] = cycler.cycler('color', color_cycle) + +matplotlib.use('Agg') + +all_wells = [row + str(i).zfill(2) for row in string.ascii_uppercase[:16] + for i in range(1, 25)] + + +def get_git_revision_hash() -> str: + #  Requires git to be installed + return subprocess.check_output(['git', 'rev-parse', 'HEAD']).decode('ascii').strip() + + +def main(): + parser = argparse.ArgumentParser() + parser.add_argument('data_directory') + parser.add_argument('-c', '--no_cpus', default=1, type=int) + parser.add_argument('--output_dir') + parser.add_argument('-w', '--wells', nargs='+') + parser.add_argument('--protocols', nargs='+') + parser.add_argument('--reversal_spread_threshold', type=float, default=10) + parser.add_argument('--export_failed', action='store_true') + parser.add_argument('--selection_file') + parser.add_argument('--figsize', nargs=2, type=int, default=[5, 8]) + parser.add_argument('--debug', action='store_true') + parser.add_argument('--log_level', default='INFO') + parser.add_argument('--Erev', default=-90.71, type=float) + parser.add_argument('--output_traces', action='store_true', + help="When true output raw and processed traces as .csv files") + + args = parser.parse_args() + + logging.basicConfig(level=args.log_level) + + if args.output_dir is None: + args.output_dir = os.path.join('output', 'hergqc') + + if not os.path.exists(args.output_dir): + os.makedirs(args.output_dir) + + with open(os.path.join(args.output_dir, 'info.txt'), 'w') as description_fout: + git_hash = get_git_revision_hash() + datetimestr = str(datetime.datetime.now()) + description_fout.write(f"Date: {datetimestr}\n") + description_fout.write(f"Commit {git_hash}\n") + command = " ".join(sys.argv) + description_fout.write(f"Command: {command}\n") + + spec = importlib.util.spec_from_file_location( + 'export_config', + os.path.join(args.data_directory, + 'export_config.py')) + + if args.wells is None: + args.wells = all_wells + wells = args.wells + + else: + wells = args.wells + + # Import and exec config file + global export_config + export_config = importlib.util.module_from_spec(spec) + + sys.modules['export_config'] = export_config + spec.loader.exec_module(export_config) + + export_config.savedir = args.output_dir + + args.saveID = export_config.saveID + args.savedir = export_config.savedir + args.D2S = export_config.D2S + args.D2SQC = export_config.D2S_QC + + protocols_regex = \ + r'^([a-z|A-Z|_|0-9| |\-|\(|\)]+)_([0-9][0-9]\.[0-9][0-9]\.[0-9][0-9])$' + + protocols_regex = re.compile(protocols_regex) + + res_dict = {} + for dirname in os.listdir(args.data_directory): + dirname = os.path.basename(dirname) + match = protocols_regex.match(dirname) + + if match is None: + continue + + protocol_name = match.group(1) + + if protocol_name not in export_config.D2S\ + and protocol_name not in export_config.D2S_QC: + continue + + # map name to new name using export_config + # savename = export_config.D2S[protocol_name] + time = match.group(2) + + if protocol_name not in res_dict: + res_dict[protocol_name] = [] + + res_dict[protocol_name].append(time) + + readnames, savenames, times_list = [], [], [] + + combined_dict = {**export_config.D2S, **export_config.D2S_QC} + + # Select QC protocols and times + for protocol in res_dict: + if protocol not in export_config.D2S_QC: + continue + + times = sorted(res_dict[protocol]) + + savename = export_config.D2S_QC[protocol] + + if len(times) == 2: + savenames.append(savename) + readnames.append(protocol) + times_list.append(times) + + elif len(times) == 4: + savenames.append(savename) + readnames.append(protocol) + times_list.append([times[0], times[2]]) + + # Make seperate savename for protocol repeat + savename = combined_dict[protocol] + '_2' + assert savename not in export_config.D2S.values() + savenames.append(savename) + times_list.append([times[1], times[3]]) + readnames.append(protocol) + + qc_dfs = [] + # Do QC which requires both repeats + passed_qc_dict = {} + for time_strs, (readname, savename) in zip(times_list, export_config.D2S_QC.items()): + for well in wells: + passed, qc_df = run_qc(readname, savename, well, time_strs, + args) + qc_dfs.append(qc_df) + passed_qc_dict[well] = True + + qc_df = pd.concat(qc_dfs, ignore_index=True) + + + # Write qc_df to file + qc_df.to_csv(os.path.join(args.savedir, 'QC-%s.csv' % args.saveID)) + + # Write data to JSON file + qc_df.to_json(os.path.join(args.savedir, 'QC-%s.json' % args.saveID), + orient='records') + + # Overwrite old files + for protocol in list(export_config.D2S_QC.values()): + fname = os.path.join(args.savedir, 'selected-%s-%s.txt' % (args.saveID, protocol)) + with open(fname, 'w') as fout: + pass + + chrono_dict = {times[0]: prot for prot, times in zip(savenames, times_list)} + + with open(os.path.join(args.output_dir, 'chrono.txt'), 'w') as fout: + for key in sorted(chrono_dict): + val = chrono_dict[key] + #  Output order of protocols + fout.write(val) + fout.write('\n') + + for time_strs, (readname, savename) in zip(times_list, export_config.D2S_QC.items()): + for well in wells: + qc_df = extract_protocol(readname, savename, well, time_strs, + args) + qc_dfs.append(qc_df) + + qc_styled_df = create_qc_table(qc_df) + logging.info(qc_styled_df) + + qc_styled_df.to_excel(os.path.join(args.output_dir, 'qc_table.xlsx')) + qc_styled_df.to_latex(os.path.join(args.output_dir, 'qc_table.tex')) + + # Save in csv format + qc_df.to_csv(os.path.join(args.savedir, 'QC-%s.csv' % args.saveID)) + + # Write data to JSON file + qc_df.to_json(os.path.join(args.savedir, 'QC-%s.json' % args.saveID), + orient='records') + + #  Load only QC vals. TODO use a new variabile name to avoid confusion + qc_vals_df = qc_df[['well', 'sweep', 'protocol', 'Rseal', 'Cm', 'Rseries']].copy() + qc_vals_df['drug'] = 'before' + qc_vals_df.to_csv(os.path.join(args.output_dir, 'qc_vals_df.csv')) + + qc_df.to_csv(os.path.join(args.output_dir, 'subtraction_qc.csv')) + + with open(os.path.join(args.output_dir, 'passed_wells.txt'), 'w') as fout: + for well, passed in passed_qc_dict.items(): + if passed: + fout.write(well) + fout.write('\n') + + +def create_qc_table(qc_df): + if len(qc_df.index) == 0: + return None + + if 'Unnamed: 0' in qc_df: + qc_df = qc_df.drop('Unnamed: 0', axis='columns') + + qc_criteria = list(qc_df.drop(['protocol', 'well'], axis='columns').columns) + + def agg_func(x): + x = x.values.flatten().astype(bool) + return bool(np.all(x)) + + qc_df[qc_criteria] = qc_df[qc_criteria].apply(lambda column: [elem == 'True' or elem is True for elem in column]) + + qc_df['protocol'] = ['staircaseramp1_2' if p == 'staircaseramp2' else p + for p in qc_df.protocol] + + print(qc_df.protocol.unique()) + + fails_dict = {} + no_wells = 384 + + dfs = [] + protocol_headings = ['staircaseramp1', 'staircaseramp1_2', 'all'] + for protocol in protocol_headings: + fails_dict = {} + for crit in sorted(qc_criteria) + ['all']: + if protocol != 'all': + sub_df = qc_df[qc_df.protocol == protocol].copy() + else: + sub_df = qc_df.copy() + + agg_dict = {crit: agg_func for crit in qc_criteria} + if crit != 'all': + col = sub_df.groupby('well').agg(agg_dict).reset_index()[crit] + vals = col.values.flatten() + n_passed = vals.sum() + else: + excluded = [crit for crit in qc_criteria + if 'all' in crit or 'spread' in crit or 'bookend' in crit] + if protocol == 'all': + excluded = [] + crit_included = [crit for crit in qc_criteria if crit not in excluded] + + col = sub_df.groupby('well').agg(agg_dict).reset_index() + n_passed = np.sum(np.all(col[crit_included].values, axis=1).flatten()) + + crit = re.sub('_', r'\_', crit) + fails_dict[crit] = (crit, no_wells - n_passed) + + new_df = pd.DataFrame.from_dict(fails_dict, orient='index', + columns=['crit', 'wells failing']) + new_df['protocol'] = protocol + new_df.set_index('crit') + dfs.append(new_df) + + ret_df = pd.concat(dfs, ignore_index=True) + + ret_df['wells failing'] = ret_df['wells failing'].astype(int) + + ret_df['protocol'] = pd.Categorical(ret_df['protocol'], + categories=protocol_headings, + ordered=True) + + return ret_df + + +def run_qc(readname, savename, well, time_strs, args): + + assert len(time_strs) == 2 + filepath_before = os.path.join(args.data_directory, + f"{readname}_{time_strs[0]}") + json_file_before = f"{readname}_{time_strs[0]}" + + filepath_after = os.path.join(args.data_directory, + f"{readname}_{time_strs[1]}") + json_file_after = f"{readname}_{time_strs[1]}" + + logging.debug(f"loading {json_file_after} and {json_file_before}") + + before_trace = SyncropatchTrace(filepath_before, + json_file_before) + + after_trace = SyncropatchTrace(filepath_after, + json_file_after) + + assert before_trace.sampling_rate == after_trace.sampling_rate + + sampling_rate = before_trace.sampling_rate + + savedir = args.output_dir + if not os.path.exists(savedir): + os.makedirs(savedir) + + before_voltage = before_trace.get_voltage() + after_voltage = after_trace.get_voltage() + + # Assert that protocols are exactly the same + assert np.all(before_voltage == after_voltage) + + voltage = before_voltage + + sweeps = [0, 1] + raw_before_all = before_trace.get_trace_sweeps(sweeps) + raw_after_all = after_trace.get_trace_sweeps(sweeps) + + selected_wells = [] + + hergqc = hERGQC(sampling_rate=sampling_rate, + # plot_dir=plot_dir, + voltage=before_voltage) + + plot_dir = os.path.join(savedir, "debug", f"debug_{well}_{savename}") + + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + qc_before = before_trace.get_onboard_QC_values() + qc_after = after_trace.get_onboard_QC_values() + + # Check if any cell first! + if (None in qc_before[well][0]) or (None in qc_after[well][0]): + no_cell = True + return False, pd.DataFrame() + + else: + no_cell = False + + nsweeps = before_trace.NofSweeps + assert after_trace.NofSweeps == nsweeps + + before_currents_corrected = np.empty((nsweeps, before_trace.NofSamples)) + after_currents_corrected = np.empty((nsweeps, after_trace.NofSamples)) + + before_currents = np.empty((nsweeps, before_trace.NofSamples)) + after_currents = np.empty((nsweeps, after_trace.NofSamples)) + + # Get ramp times from protocol description + voltage_protocol = VoltageProtocol.from_voltage_trace(voltage, + before_trace.get_times()) + + #  Find start of leak section + desc = voltage_protocol.get_all_sections() + ramp_locs = np.argwhere(desc[:, 2] != desc[:, 3]).flatten() + tstart = desc[ramp_locs[0], 0] + tend = voltage_protocol.get_ramps()[0][1] + + times = before_trace.get_times() + + ramp_bounds = [np.argmax(times > tstart), np.argmax(times > tend)] + + assert after_trace.NofSamples == before_trace.NofSamples + + for sweep in range(nsweeps): + before_raw = np.array(raw_before_all[well])[sweep, :] + after_raw = np.array(raw_after_all[well])[sweep, :] + + before_params1, before_leak = fit_linear_leak(before_raw, + voltage, + times, + *ramp_bounds, + save_fname=f"{well}-sweep{sweep}-before.png", + output_dir=savedir) + + after_params1, after_leak = fit_linear_leak(after_raw, + voltage, + times, + *ramp_bounds, + save_fname=f"{well}-sweep{sweep}-after.png", + output_dir=savedir) + + before_currents_corrected[sweep, :] = before_raw - before_leak + after_currents_corrected[sweep, :] = after_raw - after_leak + + before_currents[sweep, :] = before_raw + after_currents[sweep, :] = after_raw + + logging.info(f"{well} {savename}\n----------") + logging.info(f"sampling_rate is {sampling_rate}") + + voltage_steps = [tend + for tstart, tend, vstart, vend in + voltage_protocol.get_all_sections() if vend == vstart] + + # Run QC with raw currents + passed, QC = hergqc.run_qc(voltage_steps, times, + before_currents_corrected, + after_currents_corrected, + np.array(qc_before[well])[0, :], + np.array(qc_after[well])[0, :], nsweeps) + + QC = list(QC) + df_rows = [[well] + list(QC)] + + selected = np.all(QC) and not no_cell + if selected: + selected_wells.append(well) + + # Save subtracted current in csv file + header = "\"current\"" + + for i in range(nsweeps): + + savepath = os.path.join(savedir, + f"{args.saveID}-{savename}-{well}-sweep{i}.csv") + subtracted_current = before_currents_corrected[i, :] - after_currents_corrected[i, :] + + if args.output_traces: + if not os.path.exists(savedir): + os.makedirs(savedir) + + np.savetxt(savepath, subtracted_current, delimiter=',', + comments='', header=header) + + column_labels = ['well', 'qc1.rseal', 'qc1.cm', 'qc1.rseries', 'qc2.raw', + 'qc2.subtracted', 'qc3.raw', 'qc3.E4031', 'qc3.subtracted', + 'qc4.rseal', 'qc4.cm', 'qc4.rseries', 'qc5.staircase', + 'qc5.1.staircase', 'qc6.subtracted', 'qc6.1.subtracted', + 'qc6.2.subtracted'] + + df = pd.DataFrame(np.array(df_rows), columns=column_labels) + + missing_wells_dfs = [] + # Add onboard qc to dataframe + for well in args.wells: + if well not in df['well'].values: + onboard_qc_df = pd.DataFrame([[well] + [False for col in + list(df)[1:]]], + columns=list(df)) + missing_wells_dfs.append(onboard_qc_df) + df = pd.concat([df] + missing_wells_dfs, ignore_index=True) + + df['protocol'] = savename + + return passed, df + + +def extract_protocol(readname, savename, well, time_strs, args): + logging.info(f"extracting {savename}") + savedir = args.output_dir + saveID = args.saveID + + traces_dir = os.path.join(savedir, 'traces') + + if not os.path.exists(traces_dir): + try: + os.makedirs(traces_dir) + except FileExistsError: + pass + + row_dict = {} + + subtraction_plots_dir = os.path.join(savedir, 'subtraction_plots') + + if not os.path.isdir(subtraction_plots_dir): + try: + os.makedirs(subtraction_plots_dir) + except FileExistsError: + pass + + logging.info(f"Exporting {readname} as {savename}") + + filepath_before = os.path.join(args.data_directory, + f"{readname}_{time_strs[0]}") + filepath_after = os.path.join(args.data_directory, + f"{readname}_{time_strs[1]}") + json_file_before = f"{readname}_{time_strs[0]}" + json_file_after = f"{readname}_{time_strs[1]}" + before_trace = SyncropatchTrace(filepath_before, + json_file_before) + after_trace = SyncropatchTrace(filepath_after, + json_file_after) + + voltage_protocol = before_trace.get_voltage_protocol() + times = before_trace.get_times() + voltages = before_trace.get_voltage() + + #  Find start of leak section + desc = voltage_protocol.get_all_sections() + ramp_bounds = detect_ramp_bounds(times, desc) + tstart, tend = ramp_bounds + + nsweeps_before = before_trace.NofSweeps = 2 + nsweeps_after = after_trace.NofSweeps = 2 + + assert nsweeps_before == nsweeps_after + + # Time points + times_before = before_trace.get_times() + times_after = after_trace.get_times() + + try: + assert all(np.abs(times_before - times_after) < 1e-8) + except Exception as exc: + logging.warning(f"Exception thrown when handling {savename}: ", str(exc)) + return + + header = "\"current\"" + + qc_before = before_trace.get_onboard_QC_values() + qc_after = after_trace.get_onboard_QC_values() + qc_vals_all = before_trace.get_onboard_QC_values() + + voltage_before = before_trace.get_voltage() + voltage_after = after_trace.get_voltage() + + assert len(voltage_before) == len(voltage_after) + assert len(voltage_before) == len(times_before) + assert len(voltage_after) == len(times_after) + + voltage = voltage_before + + voltage_df = pd.DataFrame(np.vstack((times_before.flatten(), + voltage.flatten())).T, + columns=['time', 'voltage']) + + if not os.path.exists(os.path.join(traces_dir, + f"{saveID}-{savename}-voltages.csv")): + voltage_df.to_csv(os.path.join(traces_dir, + f"{saveID}-{savename}-voltages.csv")) + + np.savetxt(os.path.join(traces_dir, f"{saveID}-{savename}-times.csv"), + times_before) + + # plot subtraction + fig = plt.figure(figsize=args.figsize, layout='constrained') + + reversal_plot_dir = os.path.join(savedir, 'reversal_plots') + + rows = [] + + + before_current = before_trace.get_trace_sweeps()[well] + after_current = after_trace.get_trace_sweeps()[well] + + before_leak_currents = [] + after_leak_currents = [] + + out_dir = os.path.join(savedir, + f"{saveID}-{savename}-leak_fit-before") + + for sweep in range(before_current.shape[0]): + row_dict = { + 'well': well, + 'sweep': sweep, + 'protocol': savename + } + + qc_vals = qc_vals_all[well][sweep] + if qc_vals is None: + continue + if len(qc_vals) == 0: + continue + + row_dict['Rseal'] = qc_vals[0] + row_dict['Cm'] = qc_vals[1] + row_dict['Rseries'] = qc_vals[2] + + before_params, before_leak = fit_linear_leak(before_current[sweep, :], + voltages, times, + *ramp_bounds, + output_dir=out_dir, + save_fname=f"{well}_sweep{sweep}.png" + ) + + before_leak_currents.append(before_leak) + + out_dir = os.path.join(savedir, + f"{saveID}-{savename}-leak_fit-after") + # Convert linear regression parameters into conductance and reversal + row_dict['gleak_before'] = before_params[1] + row_dict['E_leak_before'] = -before_params[0] / before_params[1] + + after_params, after_leak = fit_linear_leak(after_current[sweep, :], + voltages, times, + *ramp_bounds, + save_fname=f"{well}_sweep{sweep}.png", + output_dir=out_dir) + + after_leak_currents.append(after_leak) + + # Convert linear regression parameters into conductance and reversal + row_dict['gleak_after'] = after_params[1] + row_dict['E_leak_after'] = -after_params[0] / after_params[1] + + subtracted_trace = before_current[sweep, :] - before_leak\ + - (after_current[sweep, :] - after_leak) + after_corrected = after_current[sweep, :] - after_leak + before_corrected = before_current[sweep, :] - before_leak + + E_rev_before = infer_reversal_potential(before_corrected, times, + desc, voltages, plot=True, + output_path=os.path.join(reversal_plot_dir, + f"{well}_{savename}_sweep{sweep}_before"), + known_Erev=args.Erev) + + E_rev_after = infer_reversal_potential(after_corrected, times, + desc, voltages, + plot=True, + output_path=os.path.join(reversal_plot_dir, + f"{well}_{savename}_sweep{sweep}_after"), + known_Erev=args.Erev) + + E_rev = infer_reversal_potential(subtracted_trace, times, desc, + voltages, plot=True, + output_path=os.path.join(reversal_plot_dir, + f"{well}_{savename}_sweep{sweep}_subtracted"), + known_Erev=args.Erev) + + row_dict['R_leftover'] =\ + np.sqrt(np.sum((after_corrected)**2)/(np.sum(before_corrected**2))) + + QC_R_leftover = np.all(row_dict['R_leftover'] < 0.5) + row_dict['QC.R_leftover'] = QC_R_leftover + + row_dict['E_rev'] = E_rev + row_dict['E_rev_before'] = E_rev_before + row_dict['E_rev_after'] = E_rev_after + + row_dict['QC.Erev'] = E_rev < -50 and E_rev > -120 + + # Check QC6 for each protocol (not just the staircase) + plot_dir = os.path.join(savedir, 'debug') + + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC(sampling_rate=before_trace.sampling_rate, + plot_dir=plot_dir, + n_sweeps=before_trace.NofSweeps) + + times = before_trace.get_times() + voltage = before_trace.get_voltage() + voltage_protocol = before_trace.get_voltage_protocol() + + voltage_steps = [tstart + for tstart, tend, vstart, vend in + voltage_protocol.get_all_sections() if vend == vstart] + + current = hergqc.filter_capacitive_spikes(before_corrected - after_corrected, + times, voltage_steps) + + if args.output_traces: + out_fname = os.path.join(traces_dir, + f"{saveID}-{savename}-{well}-sweep{sweep}-subtracted.csv") + + np.savetxt(out_fname, subtracted_trace.flatten()) + rows.append(row_dict) + + extract_df = pd.DataFrame.from_dict(rows) + logging.debug(extract_df) + + times = before_trace.get_times() + voltages = before_trace.get_voltage() + + before_current_all = before_trace.get_trace_sweeps() + after_current_all = after_trace.get_trace_sweeps() + + before_current = before_current_all[well] + after_current = after_current_all[well] + + sub_df = extract_df[extract_df.well == well] + + if len(sub_df.index) == 0: + return pd.DataFrame() + + sweeps = sorted(list(sub_df.sweep.unique())) + + do_subtraction_plot(fig, times, sweeps, before_current, after_current, + voltages, ramp_bounds, well=well, + protocol=savename) + + fig.savefig(os.path.join(subtraction_plots_dir, + f"{saveID}-{savename}-{well}-sweep{sweep}-subtraction")) + fig.clf() + plt.close(fig) + + protocol_dir = os.path.join(traces_dir, 'protocols') + if not os.path.exists(protocol_dir): + try: + os.makedirs(protocol_dir) + except FileExistsError: + pass + + # extract protocol + protocol = before_trace.get_voltage_protocol() + protocol.export_txt(os.path.join(protocol_dir, + f"{saveID}-{savename}.txt")) + + json_protocol = before_trace.get_voltage_protocol_json() + + with open(os.path.join(protocol_dir, f"{saveID}-{savename}.json"), 'w') as fout: + json.dump(json_protocol, fout) + + return extract_df + + +if __name__ == '__main__': + main() From 0e9a781d93c929580426a3965553352bbc942270 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 18 Nov 2024 16:27:48 +0000 Subject: [PATCH 15/33] #21 Add test_qc5 --- pcpostprocess/hergQC.py | 9 +++++---- tests/test_herg_qc.py | 30 ++++++++++++++++++++++++++++-- 2 files changed, 33 insertions(+), 6 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index 418a5dc..9356173 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -8,7 +8,7 @@ import scipy.stats -class QCDict(object): +class QCDict: labels = [ "qc1.rseal", @@ -68,7 +68,7 @@ def all_passed(self): return all(self.passed_list()) -class hERGQC(object): +class hERGQC: def __init__(self, sampling_rate=5, plot_dir=None, voltage=np.array([]), n_sweeps=None, removal_time=5): @@ -399,7 +399,8 @@ def qc5(self, recording1, recording2, win=None, label=''): i, f = 0, None if self.plot_dir and self._debug: - plt.axvspan(win[0], win[1], color='grey', alpha=.1) + if win is not None: + plt.axvspan(win[0], win[1], color='grey', alpha=.1) plt.plot(recording1, label='recording1') plt.plot(recording2, label='recording2') plt.savefig(os.path.join(self.plot_dir, f"qc5_{label}")) @@ -429,7 +430,7 @@ def qc5_1(self, recording1, recording2, win=None, label=''): i, f = 0, -1 if self.plot_dir and self._debug: - if win: + if win is not None: plt.axvspan(win[0], win[1], color='grey', alpha=.1) fig = plt.figure() ax = fig.subplots() diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index e5428bf..5b47896 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -144,7 +144,7 @@ def test_qc3(self): for i, expected in test_matrix: recording0 = np.asarray([0, 1] * 1000) - recording1 = np.asarray(recording0 + i) + recording1 = recording0 + i result = hergqc.qc3(recording0, recording1) self.assertEqual(result[0], expected, f"({result[1]})") @@ -181,8 +181,34 @@ def passed(result): # TODO: Test on select data def test_qc5(self): + plot_dir = os.path.join(self.output_dir, "test_qc5") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC( + sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage + ) + + # qc5 checks that the maximum current during the second half of the + # staircase changes by at least 75% of the raw trace after E-4031 addition + test_matrix = [ + (-2.0, True), + (-1.0, True), + (-0.75, True), + (-0.5, False), + (0.0, False), + (0.5, False), + (0.75, False), + (1.0, False), + ] + + for i, expected in test_matrix: + recording0 = np.asarray([0, 1] * 1000) + recording1 = recording0 + i + result = hergqc.qc5(recording0, recording1) + self.assertEqual(result[0], expected, f"({result[1]})") + # TODO: Test on select data - pass def test_qc6(self): # TODO: Test on select data From 69ba7ba9177c8c86d5aa27efba1f982275c20e18 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 18 Nov 2024 16:44:00 +0000 Subject: [PATCH 16/33] #21 Add test_qc5_1 --- pcpostprocess/hergQC.py | 2 -- tests/test_herg_qc.py | 38 ++++++++++++++++++++++++++++++++++++-- 2 files changed, 36 insertions(+), 4 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index 9356173..ceb04e3 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -410,8 +410,6 @@ def qc5(self, recording1, recording2, win=None, label=''): max_diff = recording1[i:f][wherepeak] - recording2[i:f][wherepeak] max_diffc = self.max_diffc * recording1[i:f][wherepeak] - logging.debug(f"qc5: max_diff = {max_diff}, max_diffc = {max_diffc}") - 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}") diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 5b47896..af924ec 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -146,7 +146,7 @@ def test_qc3(self): recording0 = np.asarray([0, 1] * 1000) recording1 = recording0 + i result = hergqc.qc3(recording0, recording1) - self.assertEqual(result[0], expected, f"({result[1]})") + self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data @@ -206,7 +206,41 @@ def test_qc5(self): recording0 = np.asarray([0, 1] * 1000) recording1 = recording0 + i result = hergqc.qc5(recording0, recording1) - self.assertEqual(result[0], expected, f"({result[1]})") + self.assertEqual(result[0], expected, f"({i}: {result[1]})") + + # TODO: Test on select data + + def test_qc5_1(self): + plot_dir = os.path.join(self.output_dir, "test_qc5") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC( + sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage + ) + + # qc5_1 checks that the RMSD to zero of staircase protocol changes + # by at least 50% of the raw trace after E-4031 addition. + test_matrix = [ + (0.1, True), + (0.2, True), + (0.3, True), + (0.4, True), + (0.49, True), + (0.5, True), + (0.51, False), + (0.6, False), + (0.7, False), + (0.8, False), + (0.9, False), + (1.0, False), + ] + + for i, expected in test_matrix: + recording0 = np.asarray([0, 1] * 1000) + recording1 = recording0 * i + result = hergqc.qc5_1(recording0, recording1) + self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data From b976f5f2cf9283e7f7149ada76f8097cf67a1cd8 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Mon, 18 Nov 2024 17:01:52 +0000 Subject: [PATCH 17/33] #21 Update test_qc2 --- pcpostprocess/hergQC.py | 3 ++- tests/test_herg_qc.py | 25 +++++++++++-------------- 2 files changed, 13 insertions(+), 15 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index ceb04e3..2342151 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -461,7 +461,8 @@ def qc6(self, recording1, win=None, label=''): val = np.mean(recording1[i:f]) if self.plot_dir and self._debug: - plt.axvspan(win[0], win[1], color='grey', alpha=.1) + if win is not None: + plt.axvspan(win[0], win[1], color='grey', alpha=.1) plt.plot(recording1, label='recording1') plt.savefig(os.path.join(self.plot_dir, f"qc6_{label}")) plt.clf() diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index af924ec..4561257 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -104,21 +104,17 @@ def test_qc2(self): voltage=self.voltage) # qc2 checks that raw and subtracted SNR are above a minimum threshold - recording = np.asarray([0, 1] * 500 + [0, 10] * 500) # snr = 70.75 - result = hergqc.qc2(recording) - self.assertTrue(result[0], f"({result[1]})") - - recording = np.asarray([0, 1] * 500 + [0, 6.03125] * 500) # snr = 25.02 - result = hergqc.qc2(recording) - self.assertTrue(result[0], f"({result[1]})") - - recording = np.asarray([0, 1] * 500 + [0, 6.015625] * 500) # snr = 24.88 - result = hergqc.qc2(recording) - self.assertFalse(result[0], f"({result[1]})") + test_matrix = [ + (10, True), # snr = 70.75 + (6.03125, True), # snr = 25.02 + (6.015625, False), # snr = 24.88 + (1, False), # snr = 1.0 + ] - recording = np.asarray([0, 1] * 1000) # snr = 1.0 - result = hergqc.qc2(recording) - self.assertFalse(result[0], f"({result[1]})") + for i, expected in test_matrix: + recording = np.asarray([0, 1] * 500 + [0, i] * 500) # snr = 70.75 + result = hergqc.qc2(recording) + self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data @@ -200,6 +196,7 @@ def test_qc5(self): (0.5, False), (0.75, False), (1.0, False), + (2.0, False), ] for i, expected in test_matrix: From 1b273f642855885ee40e8235ed53b80caa8415a1 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Tue, 19 Nov 2024 20:37:28 +0000 Subject: [PATCH 18/33] #21 Update qc unit tests --- pcpostprocess/hergQC.py | 26 ++--- tests/test_herg_qc.py | 211 ++++++++++++++++++++++++++++++---------- 2 files changed, 171 insertions(+), 66 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index 2342151..f5c4e4a 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -7,6 +7,7 @@ import numpy as np import scipy.stats +NOISE_LEN = 200 class QCDict: @@ -309,10 +310,10 @@ def qc2(self, recording, method=3): # Not sure if this is good... snr = scipy.stats.signaltonoise(recording) elif method == 2: - noise = np.std(recording[:200]) + noise = np.std(recording[:NOISE_LEN]) snr = (np.max(recording) - np.min(recording) - 2 * noise) / noise elif method == 3: - noise = np.std(recording[:200]) + noise = np.std(recording[:NOISE_LEN]) snr = (np.std(recording) / noise) ** 2 if snr < self.snrc or not np.isfinite(snr): @@ -328,15 +329,15 @@ def qc3(self, recording1, recording2, method=3): if method == 1: rmsdc = 2 # A/F * F elif method == 2: - noise_1 = np.std(recording1[:200]) + noise_1 = np.std(recording1[:NOISE_LEN]) peak_1 = (np.max(recording1) - noise_1) - noise_2 = np.std(recording2[:200]) + noise_2 = np.std(recording2[:NOISE_LEN]) peak_2 = (np.max(recording2) - noise_2) rmsdc = max(np.mean([peak_1, peak_2]) * 0.1, np.mean([noise_1, noise_2]) * 5) elif method == 3: - noise_1 = np.std(recording1[:200]) - noise_2 = np.std(recording2[:200]) + noise_1 = np.std(recording1[:NOISE_LEN]) + noise_2 = np.std(recording2[:NOISE_LEN]) rmsd0_1 = np.sqrt(np.mean((recording1) ** 2)) rmsd0_2 = np.sqrt(np.mean((recording2) ** 2)) rmsdc = max(np.mean([rmsd0_1, rmsd0_2]) * self.rmsd0c, @@ -348,7 +349,7 @@ def qc3(self, recording1, recording2, method=3): else: result = True - return (result, rmsd) + return (result, rmsdc - rmsd) def qc4(self, rseals, cms, rseriess): # Check R_seal, C_m, R_series stability @@ -417,7 +418,7 @@ def qc5(self, recording1, recording2, win=None, label=''): else: result = True - return (result, max_diff) + return (result, max_diffc - max_diff) def qc5_1(self, recording1, recording2, win=None, label=''): # Check RMSD_0 drops after E-4031 application @@ -450,7 +451,7 @@ def qc5_1(self, recording1, recording2, win=None, label=''): else: result = True - return (result, rmsd0_diff) + return (result, rmsd0_diffc - rmsd0_diff) def qc6(self, recording1, win=None, label=''): # Check subtracted staircase +40mV step up is non negative @@ -458,7 +459,6 @@ def qc6(self, recording1, win=None, label=''): i, f = win else: i, f = 0, -1 - val = np.mean(recording1[i:f]) if self.plot_dir and self._debug: if win is not None: @@ -467,8 +467,10 @@ def qc6(self, recording1, win=None, label=''): plt.savefig(os.path.join(self.plot_dir, f"qc6_{label}")) plt.clf() + val = np.mean(recording1[i:f]) + # valc = -0.005 * np.abs(np.sqrt(np.mean((recording1) ** 2))) # or just 0 - valc = self.negative_tolc * np.std(recording1[:200]) + valc = self.negative_tolc * np.std(recording1[:NOISE_LEN]) if (val < valc) or not (np.isfinite(val) and np.isfinite(valc)): self.logger.debug(f"qc6_{label} val: {val}, valc: {valc}") @@ -476,7 +478,7 @@ def qc6(self, recording1, win=None, label=''): else: result = True - return (result, val) + return (result, valc - 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 4561257..aed37d5 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -6,7 +6,7 @@ import numpy as np from syncropatch_export.trace import Trace -from pcpostprocess.hergQC import hERGQC +from pcpostprocess.hergQC import hERGQC, NOISE_LEN class TestHergQC(unittest.TestCase): @@ -54,15 +54,14 @@ def passed(result): # qc1 checks that rseal, cm, rseries are within range rseal_lo, rseal_hi = 1e8, 1e12 rseal_mid = (rseal_lo + rseal_hi) / 2 + rseal_tol = 0.1 cm_lo, cm_hi = 1e-12, 1e-10 cm_mid = (cm_lo + cm_hi) / 2 + cm_tol = 1e-13 rseries_lo, rseries_hi = 1e6, 2.5e7 rseries_mid = (rseries_lo + rseries_hi) / 2 - - rseal_tol = 0.1 - cm_tol = 1e-13 rseries_tol = 0.1 test_matrix = [ @@ -105,14 +104,17 @@ def test_qc2(self): # qc2 checks that raw and subtracted SNR are above a minimum threshold test_matrix = [ - (10, True), # snr = 70.75 - (6.03125, True), # snr = 25.02 - (6.015625, False), # snr = 24.88 - (1, False), # snr = 1.0 + (40, True), # snr = 130286 + (10, True), # snr = 8082 + (1, True), # snr = 74 + (0.601, True), # snr = 25.07 + (0.6, False), # snr = 24.98 + (0.5, False), # snr = 17 + (0.1, False), # snr = 0.5 ] for i, expected in test_matrix: - recording = np.asarray([0, 1] * 500 + [0, i] * 500) # snr = 70.75 + recording = np.asarray([0, 0.1] * (NOISE_LEN//2) + [i] * 500) result = hergqc.qc2(recording) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -128,20 +130,37 @@ def test_qc3(self): voltage=self.voltage) # qc3 checks that rmsd of two sweeps are similar + + # Test with same noise, different signal test_matrix = [ - (0, True), - (1, True), - (2, True), - (3, True), - (4, False), - (5, False), - (6, False), + (-10, False), + (-9, False), + (-8, False), # rmsdc - rmsd = -0.6761186497920804 + (-7, True), # rmsdc - rmsd = 0.25355095037585507 + (0, True), # rmsdc - rmsd = 6.761238263598085 + (8, True), # rmsdc - rmsd = 0.6761272774054383 + (9, False), # rmsdc - rmsd = -0.08451158778363332 + (10, False), ] + recording1 = np.asarray([0, 0.1] * (NOISE_LEN//2) + [40] * 500) for i, expected in test_matrix: - recording0 = np.asarray([0, 1] * 1000) - recording1 = recording0 + i - result = hergqc.qc3(recording0, recording1) + recording2 = np.asarray([0, 0.1] * (NOISE_LEN//2) + [40 + i] * 500) + result = hergqc.qc3(recording1, recording2) + self.assertEqual(result[0], expected, f"({i}: {result[1]})") + + # Test with same signal, different noise + # TODO: Find failing example + test_matrix = [ + (10, True), + (100, True), + (1000, True), + ] + + recording1 = np.asarray([0, 0.1] * (NOISE_LEN//2) + [40] * 500) + for i, expected in test_matrix: + recording2 = np.asarray([0, 0.1 * i] * (NOISE_LEN//2) + [40] * 500) + result = hergqc.qc3(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data @@ -159,19 +178,85 @@ def passed(result): ) # qc4 checks that rseal, cm, rseries are similar before/after E-4031 change + r_lo, r_hi = 1e6, 3e7 + c_lo, c_hi = 1e-12, 1e-10 + + # Test rseals + cms = [c_lo, c_lo] + rseriess = [r_lo, r_lo] + + test_matrix = [ + (1.1, True), + (3.0, True), + (3.5, False), + (5.0, False), + ] + + for i, expected in test_matrix: + rseals = [r_lo, i * r_lo] + self.assertEqual( + passed(hergqc.qc4(rseals, cms, rseriess)), + expected, + f"({i}: {rseals}, {cms}, {rseriess})", + ) + + rseals = [r_hi, i * r_hi] + self.assertEqual( + passed(hergqc.qc4(rseals, cms, rseriess)), + expected, + f"({i}: {rseals}, {cms}, {rseriess})", + ) + + # Test cms + rseals = [r_lo, r_lo] + rseriess = [r_lo, r_lo] + test_matrix = [ - [([1, 1], [1, 1], [1, 1]), True], - [([1, 2], [1, 2], [1, 2]), True], - [([1, 3], [1, 3], [1, 3]), True], - [([1, 4], [1, 4], [1, 4]), False], - [([1, 5], [1, 5], [1, 5]), False], + (1.1, True), + (3.0, True), + (3.5, False), + (5.0, False), ] - for (rseals, cms, rseriess), expected in test_matrix: + for i, expected in test_matrix: + cms = [c_lo, i * c_lo] self.assertEqual( passed(hergqc.qc4(rseals, cms, rseriess)), expected, - f"({rseals}, {cms}, {rseriess})", + f"({i}: {rseals}, {cms}, {rseriess})", + ) + + cms = [c_hi, i * c_hi] + self.assertEqual( + passed(hergqc.qc4(rseals, cms, rseriess)), + expected, + f"({i}: {rseals}, {cms}, {rseriess})", + ) + + # Test rseriess + cms = [c_lo, c_lo] + rseals = [r_lo, r_lo] + + test_matrix = [ + (1.1, True), + (3.0, True), + (3.5, False), + (5.0, False), + ] + + for i, expected in test_matrix: + rseriess = [r_lo, i * r_lo] + self.assertEqual( + passed(hergqc.qc4(rseals, cms, rseriess)), + expected, + f"({i}: {rseals}, {cms}, {rseriess})", + ) + + rseriess = [r_hi, i * r_hi] + self.assertEqual( + passed(hergqc.qc4(rseals, cms, rseriess)), + expected, + f"({i}: {rseals}, {cms}, {rseriess})", ) # TODO: Test on select data @@ -185,30 +270,28 @@ def test_qc5(self): sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage ) - # qc5 checks that the maximum current during the second half of the + # qc5 checks that the maximum current during the second half of the # staircase changes by at least 75% of the raw trace after E-4031 addition test_matrix = [ - (-2.0, True), (-1.0, True), - (-0.75, True), - (-0.5, False), - (0.0, False), + (0.1, True), + (0.2, True), # max_diffc - max_diff = -0.5 + (0.25, True), # max_diffc - max_diff = 0.0 + (0.3, False), # max_diffc - max_diff = 0.5 (0.5, False), - (0.75, False), (1.0, False), - (2.0, False), ] + recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) for i, expected in test_matrix: - recording0 = np.asarray([0, 1] * 1000) - recording1 = recording0 + i - result = hergqc.qc5(recording0, recording1) + recording2 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) + result = hergqc.qc5(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data def test_qc5_1(self): - plot_dir = os.path.join(self.output_dir, "test_qc5") + plot_dir = os.path.join(self.output_dir, "test_qc5_1") if not os.path.exists(plot_dir): os.makedirs(plot_dir) @@ -216,34 +299,54 @@ def test_qc5_1(self): sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage ) - # qc5_1 checks that the RMSD to zero of staircase protocol changes + # qc5_1 checks that the RMSD to zero of staircase protocol changes # by at least 50% of the raw trace after E-4031 addition. test_matrix = [ - (0.1, True), - (0.2, True), - (0.3, True), - (0.4, True), - (0.49, True), - (0.5, True), - (0.51, False), - (0.6, False), - (0.7, False), - (0.8, False), - (0.9, False), - (1.0, False), + (-1.0, False), # rmsd0_diffc - rmsd0_diff = 4.2 + (-0.5, False), # rmsd0_diffc - rmsd0_diff = 0.0001 + (-0.4, True), # rmsd0_diffc - rmsd0_diff = -0.8 + (-0.1, True), # rmsd0_diffc - rmsd0_diff = -3.4 + (0.1, True), # rmsd0_diffc - rmsd0_diff = -3.4 + (0.4, True), # rmsd0_diffc - rmsd0_diff = -0.8 + (0.5, False), # rmsd0_diffc - rmsd0_diff = 0.0001 + (1.0, False), # rmsd0_diffc - rmsd0_diff = 4.2 ] + recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) for i, expected in test_matrix: - recording0 = np.asarray([0, 1] * 1000) - recording1 = recording0 * i - result = hergqc.qc5_1(recording0, recording1) + recording2 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) + result = hergqc.qc5_1(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data def test_qc6(self): + plot_dir = os.path.join(self.output_dir, "test_qc2") + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) + + hergqc = hERGQC( + sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage + ) + + # qc6 checks that the first step up to +40 mV, before the staircase, in + # the subtracted trace is bigger than -2 x estimated noise level. + test_matrix = [ + (-1, False), # valc - val = 9.9 + (-0.1, False), # valc - val = 0.9 + (-0.02, False), # valc - val = 0.1 + (-0.01, True), # valc - val = -1.38 + (-0.001, True), # valc - val = -0.09 + (0.001, True), # valc - val = -0.11 + (0.1, True), # valc - val = -1.1 + (1, True), # valc - val = -10.1 + ] + for i, expected in test_matrix: + recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) + result = hergqc.qc6(recording, win=[NOISE_LEN, -1]) + self.assertEqual(result[0], expected, f"({i}: {result[1]})") + # TODO: Test on select data - pass def test_run_qc(self): self.assertTrue(np.all(np.isfinite(self.voltage))) From de5602d43606bf3c5a5456d0e96f711458260ee1 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Tue, 19 Nov 2024 20:56:05 +0000 Subject: [PATCH 19/33] #21 Update test_qc6 --- tests/test_herg_qc.py | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index aed37d5..88799a1 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -332,17 +332,17 @@ def test_qc6(self): # qc6 checks that the first step up to +40 mV, before the staircase, in # the subtracted trace is bigger than -2 x estimated noise level. test_matrix = [ - (-1, False), # valc - val = 9.9 - (-0.1, False), # valc - val = 0.9 - (-0.02, False), # valc - val = 0.1 - (-0.01, True), # valc - val = -1.38 - (-0.001, True), # valc - val = -0.09 - (0.001, True), # valc - val = -0.11 - (0.1, True), # valc - val = -1.1 - (1, True), # valc - val = -10.1 + (-100, False), # valc - val = 9.9 + (-10, False), # valc - val = 0.9 + (-2, False), # valc - val = 0.1 + (-1, True), # valc - val = 0 + (1, True), # valc - val = -0.2 + (2, True), # valc - val = -0.3 + (10, True), # valc - val = -1.1 + (100, True), # valc - val = -10.1 ] for i, expected in test_matrix: - recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) + recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [0.1 * i] * 500) result = hergqc.qc6(recording, win=[NOISE_LEN, -1]) self.assertEqual(result[0], expected, f"({i}: {result[1]})") From cca851bde111ec4734f0a774443108540c49644e Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Wed, 20 Nov 2024 17:15:24 +0000 Subject: [PATCH 20/33] #21 test_qc1 on test data --- pcpostprocess/hergQC.py | 18 ++-- tests/test_herg_qc.py | 230 +++++++++++++++++++++------------------- 2 files changed, 135 insertions(+), 113 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index f5c4e4a..b17fd2d 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -75,8 +75,7 @@ def __init__(self, sampling_rate=5, plot_dir=None, voltage=np.array([]), n_sweeps=None, removal_time=5): # TODO docstring - if plot_dir is not None: - self.plot_dir = plot_dir + self._plot_dir = plot_dir self._n_qc = 16 @@ -134,8 +133,15 @@ def __init__(self, sampling_rate=5, plot_dir=None, voltage=np.array([]), self._debug = True - def set_trace(self, before, after, qc_vals_before, - qc_vals_after, n_sweeps): + @property + def plot_dir(self): + return self._plot_dir + + @plot_dir.setter + def plot_dir(self, path): + self._plot_dir = path + + def set_trace(self, before, after, qc_vals_before, qc_vals_after, n_sweeps): self._before = before self._qc_vals_before = qc_vals_before self._after = after @@ -154,7 +160,7 @@ def run_qc(self, voltage_steps, times, @param before is the before-drug current trace @param after is the post-drug current trace @param qc_vals_before is an array of values for each pre-drug sweep where each row is (Rseal, Cm, Rseries) - @param qc_vals_before is an array of values for each post-drug sweep where each row is (Rseal, Cm, Rseries) + @param qc_vals_after is an array of values for each post-drug sweep where each row is (Rseal, Cm, Rseries) @n_sweeps is the number of sweeps to be considered """ @@ -241,7 +247,7 @@ def run_qc(self, voltage_steps, times, QC['qc6.1.subtracted'].append(qc6_1) QC['qc6.2.subtracted'].append(qc6_2) - if self._debug: + if self.plot_dir and self._debug: fig = plt.figure(figsize=(8, 5)) ax = fig.subplots() ax.plot(times, (before - after).T, label='subtracted') diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 88799a1..dec3876 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -1,5 +1,6 @@ import logging import os +import copy import string import unittest @@ -12,44 +13,60 @@ class TestHergQC(unittest.TestCase): def setUp(self): - filepath = os.path.join('tests', 'test_data', '13112023_MW2_FF', - 'staircaseramp (2)_2kHz_15.01.07') + base_path = os.path.join("tests", "test_data", "13112023_MW2_FF") - self.all_wells = [ - lab + str(i).zfill(2) for lab in string.ascii_uppercase[:16] - for i in range(1, 25)] + label_before = "staircaseramp (2)_2kHz_15.01.07" + label_after = "staircaseramp (2)_2kHz_15.11.33" - filepath2 = os.path.join('tests', 'test_data', '13112023_MW2_FF', - 'staircaseramp (2)_2kHz_15.11.33') + path_before = os.path.join(base_path, label_before) + path_after = os.path.join(base_path, label_after) - json_file = "staircaseramp (2)_2kHz_15.01.07" - json_file2 = "staircaseramp (2)_2kHz_15.11.33" + trace_before = Trace(path_before, json_file=label_before) + trace_after = Trace(path_after, json_file=label_after) - self.output_dir = os.path.join('test_output', 'test_herg_qc') + sweeps = [0, 1] + self.trace_sweeps_before = trace_before.get_trace_sweeps(sweeps) + self.trace_sweeps_after = trace_after.get_trace_sweeps(sweeps) + + self.qc_vals_before = trace_before.get_onboard_QC_values(sweeps=sweeps) + self.qc_vals_after = trace_after.get_onboard_QC_values(sweeps=sweeps) + + self.voltage = trace_before.get_voltage() + self.times = trace_after.get_times() - if not os.path.exists(self.output_dir): - os.makedirs(self.output_dir) + sampling_rate = int(1.0 / (self.times[1] - self.times[0])) # in kHz + + #  Assume that there are no discontinuities at the start or end of ramps + voltage_protocol = trace_before.get_voltage_protocol() + self.voltage_steps = [ + tstart + for tstart, _, vstart, vend in voltage_protocol.get_all_sections() + if vend == vstart + ] - self.test_trace_before = Trace(filepath, json_file) - self.test_trace_after = Trace(filepath2, json_file2) + plot_dir = "test_output" + if not os.path.exists(plot_dir): + os.makedirs(plot_dir) - self.voltage = self.test_trace_before.get_voltage() - self.times = self.test_trace_after.get_times() + self.hergqc = hERGQC( + sampling_rate=sampling_rate, + plot_dir=plot_dir, + voltage=self.voltage, + ) - # Calculate sampling rate in (use kHz) - self.sampling_rate = int(1.0 / (self.times[1] - self.times[0])) + def test_qc_inputs(self): + self.assertTrue(np.all(np.isfinite(self.voltage))) + self.assertTrue(np.all(np.isfinite(self.times))) def test_qc1(self): def passed(result): return all([x for x, _ in result]) - plot_dir = os.path.join(self.output_dir, "test_qc1") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc1") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC(sampling_rate=self.sampling_rate, - plot_dir=plot_dir, - voltage=self.voltage) + hergqc.plot_dir = plot_dir # qc1 checks that rseal, cm, rseries are within range rseal_lo, rseal_hi = 1e8, 1e12 @@ -88,46 +105,62 @@ def passed(result): self.assertEqual( passed(hergqc.qc1(rseal, cm, rseries)), expected, - f"({rseal}, {cm}, {rseries})", + f"QC1: {rseal}, {cm}, {rseries}", ) - # TODO: Test on select data + # Test on select data + test_matrix = [ + ("A01", True), + ("A02", True), + ("A03", True), + ("A04", True), + ("A05", True), + ("A10", False), + ("A12", False), + ("A13", False), + ("A16", False), + ("A19", False), + ] + + for well, expected in test_matrix: + with self.subTest(well): + qc_vals_before = np.array(self.qc_vals_before[well])[0, :] + self.assertEqual( + passed(hergqc.qc1(*qc_vals_before)), + expected, + f"QC1: {well} (before) {qc_vals_before}", + ) def test_qc2(self): - plot_dir = os.path.join(self.output_dir, "test_qc2") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc2") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC(sampling_rate=self.sampling_rate, - plot_dir=plot_dir, - voltage=self.voltage) + hergqc.plot_dir = plot_dir # qc2 checks that raw and subtracted SNR are above a minimum threshold test_matrix = [ - (40, True), # snr = 130286 (10, True), # snr = 8082 (1, True), # snr = 74 - (0.601, True), # snr = 25.07 + (0.601, True), # snr = 25.07 (0.6, False), # snr = 24.98 (0.5, False), # snr = 17 (0.1, False), # snr = 0.5 ] for i, expected in test_matrix: - recording = np.asarray([0, 0.1] * (NOISE_LEN//2) + [i] * 500) + recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [i] * 500) result = hergqc.qc2(recording) self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data def test_qc3(self): - plot_dir = os.path.join(self.output_dir, "test_qc3") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc3") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC(sampling_rate=self.sampling_rate, - plot_dir=plot_dir, - voltage=self.voltage) + hergqc.plot_dir = plot_dir # qc3 checks that rmsd of two sweeps are similar @@ -135,17 +168,17 @@ def test_qc3(self): test_matrix = [ (-10, False), (-9, False), - (-8, False), # rmsdc - rmsd = -0.6761186497920804 - (-7, True), # rmsdc - rmsd = 0.25355095037585507 - (0, True), # rmsdc - rmsd = 6.761238263598085 - (8, True), # rmsdc - rmsd = 0.6761272774054383 - (9, False), # rmsdc - rmsd = -0.08451158778363332 + (-8, False), # rmsdc - rmsd = -0.6761186497920804 + (-7, True), # rmsdc - rmsd = 0.25355095037585507 + (0, True), # rmsdc - rmsd = 6.761238263598085 + (8, True), # rmsdc - rmsd = 0.6761272774054383 + (9, False), # rmsdc - rmsd = -0.08451158778363332 (10, False), ] - recording1 = np.asarray([0, 0.1] * (NOISE_LEN//2) + [40] * 500) + recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) for i, expected in test_matrix: - recording2 = np.asarray([0, 0.1] * (NOISE_LEN//2) + [40 + i] * 500) + recording2 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40 + i] * 500) result = hergqc.qc3(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -157,9 +190,9 @@ def test_qc3(self): (1000, True), ] - recording1 = np.asarray([0, 0.1] * (NOISE_LEN//2) + [40] * 500) + recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) for i, expected in test_matrix: - recording2 = np.asarray([0, 0.1 * i] * (NOISE_LEN//2) + [40] * 500) + recording2 = np.asarray([0, 0.1 * i] * (NOISE_LEN // 2) + [40] * 500) result = hergqc.qc3(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -169,13 +202,11 @@ def test_qc4(self): def passed(result): return all([x for x, _ in result]) - plot_dir = os.path.join(self.output_dir, "test_qc1") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc4") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC( - sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage - ) + hergqc.plot_dir = plot_dir # qc4 checks that rseal, cm, rseries are similar before/after E-4031 change r_lo, r_hi = 1e6, 3e7 @@ -262,13 +293,11 @@ def passed(result): # TODO: Test on select data def test_qc5(self): - plot_dir = os.path.join(self.output_dir, "test_qc5") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc5") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC( - sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage - ) + hergqc.plot_dir = plot_dir # qc5 checks that the maximum current during the second half of the # staircase changes by at least 75% of the raw trace after E-4031 addition @@ -291,13 +320,11 @@ def test_qc5(self): # TODO: Test on select data def test_qc5_1(self): - plot_dir = os.path.join(self.output_dir, "test_qc5_1") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc5_1") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC( - sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage - ) + hergqc.plot_dir = plot_dir # qc5_1 checks that the RMSD to zero of staircase protocol changes # by at least 50% of the raw trace after E-4031 addition. @@ -321,13 +348,11 @@ def test_qc5_1(self): # TODO: Test on select data def test_qc6(self): - plot_dir = os.path.join(self.output_dir, "test_qc2") + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_qc6") if not os.path.exists(plot_dir): os.makedirs(plot_dir) - - hergqc = hERGQC( - sampling_rate=self.sampling_rate, plot_dir=plot_dir, voltage=self.voltage - ) + hergqc.plot_dir = plot_dir # qc6 checks that the first step up to +40 mV, before the staircase, in # the subtracted trace is bigger than -2 x estimated noise level. @@ -349,56 +374,47 @@ def test_qc6(self): # TODO: Test on select data def test_run_qc(self): - self.assertTrue(np.all(np.isfinite(self.voltage))) - self.assertTrue(np.all(np.isfinite(self.times))) - - plot_dir = os.path.join(self.output_dir, - 'test_run_qc') + # Spot check a few wells; could check all, but it's time consuming. + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, "test_run_qc") if not os.path.exists(plot_dir): os.makedirs(plot_dir) + hergqc.plot_dir = plot_dir - hergqc = hERGQC(sampling_rate=self.sampling_rate, - plot_dir=plot_dir, - voltage=self.voltage) - - sweeps = [0, 1] - before = self.test_trace_before.get_trace_sweeps(sweeps) - after = self.test_trace_after.get_trace_sweeps(sweeps) - qc_vals_before = self.test_trace_before.get_onboard_QC_values(sweeps=sweeps) - qc_vals_after = self.test_trace_after.get_onboard_QC_values(sweeps=sweeps) - - # Spot check a few wells - # We could check all of the wells but it's time consuming - - test_wells = ['A01', 'A02', 'A03'] - - voltage_protocol = self.test_trace_before.get_voltage_protocol() + test_matrix = [ + ("A01", True), + ("A02", True), + ("A03", True), + ("A04", False), + ("A05", False), + ("D01", False), + ] - for well in test_wells: + for well, expected in test_matrix: with self.subTest(well): # Take values from the first sweep only - qc_vals_before_well = np.array(qc_vals_before[well])[0, :] - qc_vals_after_well = np.array(qc_vals_after[well])[0, :] - - before_well = np.array(before[well]) - after_well = np.array(after[well]) - - #  Assume that there are no discontinuities at the start or end of ramps - voltage_steps = [tstart - for tstart, tend, vstart, vend in - voltage_protocol.get_all_sections() if vend == vstart] - - QC = hergqc.run_qc(voltage_steps, - self.times, before_well, after_well, - qc_vals_before_well, - qc_vals_after_well, n_sweeps=2) + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) + + qc_vals_before = np.array(self.qc_vals_before[well])[0, :] + qc_vals_after = np.array(self.qc_vals_after[well])[0, :] + + QC = hergqc.run_qc( + voltage_steps=self.voltage_steps, + times=self.times, + before=before, + after=after, + qc_vals_before=qc_vals_before, + qc_vals_after=qc_vals_after, + n_sweeps=2, + ) logging.debug(well, QC.all_passed()) trace = "" for label, result in QC.items(): if not QC.qc_passed(label): - trace += f"{label}: {result}\n" - print(f"Testing Well, {well}") - self.assertTrue(QC.all_passed(), trace) + trace += f"{well} {label}: {result}\n" + + self.assertEqual(QC.all_passed(), expected, trace) From d615b8c5a6c1fafcd26ca3eb4f66e9686de83131 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Wed, 27 Nov 2024 17:26:22 +0000 Subject: [PATCH 21/33] #21 Remove ghost export_staircase_data from branch --- .../scripts/export_staircase_data.py | 730 ------------------ 1 file changed, 730 deletions(-) delete mode 100644 pcpostprocess/scripts/export_staircase_data.py diff --git a/pcpostprocess/scripts/export_staircase_data.py b/pcpostprocess/scripts/export_staircase_data.py deleted file mode 100644 index ddf7b20..0000000 --- a/pcpostprocess/scripts/export_staircase_data.py +++ /dev/null @@ -1,730 +0,0 @@ -import argparse -import datetime -import importlib.util -import json -import logging -import multiprocessing -import os -import string -import subprocess -import sys - -import cycler -import matplotlib -import matplotlib.pyplot as plt -import numpy as np -import pandas as pd -import regex as re -import scipy -from syncropatch_export.trace import Trace as SyncropatchTrace -from syncropatch_export.voltage_protocols import VoltageProtocol - -from pcpostprocess.detect_ramp_bounds import detect_ramp_bounds -from pcpostprocess.hergQC import hERGQC -from pcpostprocess.infer_reversal import infer_reversal_potential -from pcpostprocess.leak_correct import fit_linear_leak, get_leak_corrected -from pcpostprocess.subtraction_plots import do_subtraction_plot - -pool_kws = {'maxtasksperchild': 1} - -color_cycle = ["#5790fc", "#f89c20", "#e42536", "#964a8b", "#9c9ca1", "#7a21dd"] -plt.rcParams['axes.prop_cycle'] = cycler.cycler('color', color_cycle) - -matplotlib.use('Agg') - -all_wells = [row + str(i).zfill(2) for row in string.ascii_uppercase[:16] - for i in range(1, 25)] - - -def get_git_revision_hash() -> str: - #  Requires git to be installed - return subprocess.check_output(['git', 'rev-parse', 'HEAD']).decode('ascii').strip() - - -def main(): - parser = argparse.ArgumentParser() - parser.add_argument('data_directory') - parser.add_argument('-c', '--no_cpus', default=1, type=int) - parser.add_argument('--output_dir') - parser.add_argument('-w', '--wells', nargs='+') - parser.add_argument('--protocols', nargs='+') - parser.add_argument('--reversal_spread_threshold', type=float, default=10) - parser.add_argument('--export_failed', action='store_true') - parser.add_argument('--selection_file') - parser.add_argument('--figsize', nargs=2, type=int, default=[5, 8]) - parser.add_argument('--debug', action='store_true') - parser.add_argument('--log_level', default='INFO') - parser.add_argument('--Erev', default=-90.71, type=float) - parser.add_argument('--output_traces', action='store_true', - help="When true output raw and processed traces as .csv files") - - args = parser.parse_args() - - logging.basicConfig(level=args.log_level) - - if args.output_dir is None: - args.output_dir = os.path.join('output', 'hergqc') - - if not os.path.exists(args.output_dir): - os.makedirs(args.output_dir) - - with open(os.path.join(args.output_dir, 'info.txt'), 'w') as description_fout: - git_hash = get_git_revision_hash() - datetimestr = str(datetime.datetime.now()) - description_fout.write(f"Date: {datetimestr}\n") - description_fout.write(f"Commit {git_hash}\n") - command = " ".join(sys.argv) - description_fout.write(f"Command: {command}\n") - - spec = importlib.util.spec_from_file_location( - 'export_config', - os.path.join(args.data_directory, - 'export_config.py')) - - if args.wells is None: - args.wells = all_wells - wells = args.wells - - else: - wells = args.wells - - # Import and exec config file - global export_config - export_config = importlib.util.module_from_spec(spec) - - sys.modules['export_config'] = export_config - spec.loader.exec_module(export_config) - - export_config.savedir = args.output_dir - - args.saveID = export_config.saveID - args.savedir = export_config.savedir - args.D2S = export_config.D2S - args.D2SQC = export_config.D2S_QC - - protocols_regex = \ - r'^([a-z|A-Z|_|0-9| |\-|\(|\)]+)_([0-9][0-9]\.[0-9][0-9]\.[0-9][0-9])$' - - protocols_regex = re.compile(protocols_regex) - - res_dict = {} - for dirname in os.listdir(args.data_directory): - dirname = os.path.basename(dirname) - match = protocols_regex.match(dirname) - - if match is None: - continue - - protocol_name = match.group(1) - - if protocol_name not in export_config.D2S\ - and protocol_name not in export_config.D2S_QC: - continue - - # map name to new name using export_config - # savename = export_config.D2S[protocol_name] - time = match.group(2) - - if protocol_name not in res_dict: - res_dict[protocol_name] = [] - - res_dict[protocol_name].append(time) - - readnames, savenames, times_list = [], [], [] - - combined_dict = {**export_config.D2S, **export_config.D2S_QC} - - # Select QC protocols and times - for protocol in res_dict: - if protocol not in export_config.D2S_QC: - continue - - times = sorted(res_dict[protocol]) - - savename = export_config.D2S_QC[protocol] - - if len(times) == 2: - savenames.append(savename) - readnames.append(protocol) - times_list.append(times) - - elif len(times) == 4: - savenames.append(savename) - readnames.append(protocol) - times_list.append([times[0], times[2]]) - - # Make seperate savename for protocol repeat - savename = combined_dict[protocol] + '_2' - assert savename not in export_config.D2S.values() - savenames.append(savename) - times_list.append([times[1], times[3]]) - readnames.append(protocol) - - qc_dfs = [] - # Do QC which requires both repeats - passed_qc_dict = {} - for time_strs, (readname, savename) in zip(times_list, export_config.D2S_QC.items()): - for well in wells: - passed, qc_df = run_qc(readname, savename, well, time_strs, - args) - qc_dfs.append(qc_df) - passed_qc_dict[well] = True - - qc_df = pd.concat(qc_dfs, ignore_index=True) - - - # Write qc_df to file - qc_df.to_csv(os.path.join(args.savedir, 'QC-%s.csv' % args.saveID)) - - # Write data to JSON file - qc_df.to_json(os.path.join(args.savedir, 'QC-%s.json' % args.saveID), - orient='records') - - # Overwrite old files - for protocol in list(export_config.D2S_QC.values()): - fname = os.path.join(args.savedir, 'selected-%s-%s.txt' % (args.saveID, protocol)) - with open(fname, 'w') as fout: - pass - - chrono_dict = {times[0]: prot for prot, times in zip(savenames, times_list)} - - with open(os.path.join(args.output_dir, 'chrono.txt'), 'w') as fout: - for key in sorted(chrono_dict): - val = chrono_dict[key] - #  Output order of protocols - fout.write(val) - fout.write('\n') - - for time_strs, (readname, savename) in zip(times_list, export_config.D2S_QC.items()): - for well in wells: - qc_df = extract_protocol(readname, savename, well, time_strs, - args) - qc_dfs.append(qc_df) - - qc_styled_df = create_qc_table(qc_df) - logging.info(qc_styled_df) - - qc_styled_df.to_excel(os.path.join(args.output_dir, 'qc_table.xlsx')) - qc_styled_df.to_latex(os.path.join(args.output_dir, 'qc_table.tex')) - - # Save in csv format - qc_df.to_csv(os.path.join(args.savedir, 'QC-%s.csv' % args.saveID)) - - # Write data to JSON file - qc_df.to_json(os.path.join(args.savedir, 'QC-%s.json' % args.saveID), - orient='records') - - #  Load only QC vals. TODO use a new variabile name to avoid confusion - qc_vals_df = qc_df[['well', 'sweep', 'protocol', 'Rseal', 'Cm', 'Rseries']].copy() - qc_vals_df['drug'] = 'before' - qc_vals_df.to_csv(os.path.join(args.output_dir, 'qc_vals_df.csv')) - - qc_df.to_csv(os.path.join(args.output_dir, 'subtraction_qc.csv')) - - with open(os.path.join(args.output_dir, 'passed_wells.txt'), 'w') as fout: - for well, passed in passed_qc_dict.items(): - if passed: - fout.write(well) - fout.write('\n') - - -def create_qc_table(qc_df): - if len(qc_df.index) == 0: - return None - - if 'Unnamed: 0' in qc_df: - qc_df = qc_df.drop('Unnamed: 0', axis='columns') - - qc_criteria = list(qc_df.drop(['protocol', 'well'], axis='columns').columns) - - def agg_func(x): - x = x.values.flatten().astype(bool) - return bool(np.all(x)) - - qc_df[qc_criteria] = qc_df[qc_criteria].apply(lambda column: [elem == 'True' or elem is True for elem in column]) - - qc_df['protocol'] = ['staircaseramp1_2' if p == 'staircaseramp2' else p - for p in qc_df.protocol] - - print(qc_df.protocol.unique()) - - fails_dict = {} - no_wells = 384 - - dfs = [] - protocol_headings = ['staircaseramp1', 'staircaseramp1_2', 'all'] - for protocol in protocol_headings: - fails_dict = {} - for crit in sorted(qc_criteria) + ['all']: - if protocol != 'all': - sub_df = qc_df[qc_df.protocol == protocol].copy() - else: - sub_df = qc_df.copy() - - agg_dict = {crit: agg_func for crit in qc_criteria} - if crit != 'all': - col = sub_df.groupby('well').agg(agg_dict).reset_index()[crit] - vals = col.values.flatten() - n_passed = vals.sum() - else: - excluded = [crit for crit in qc_criteria - if 'all' in crit or 'spread' in crit or 'bookend' in crit] - if protocol == 'all': - excluded = [] - crit_included = [crit for crit in qc_criteria if crit not in excluded] - - col = sub_df.groupby('well').agg(agg_dict).reset_index() - n_passed = np.sum(np.all(col[crit_included].values, axis=1).flatten()) - - crit = re.sub('_', r'\_', crit) - fails_dict[crit] = (crit, no_wells - n_passed) - - new_df = pd.DataFrame.from_dict(fails_dict, orient='index', - columns=['crit', 'wells failing']) - new_df['protocol'] = protocol - new_df.set_index('crit') - dfs.append(new_df) - - ret_df = pd.concat(dfs, ignore_index=True) - - ret_df['wells failing'] = ret_df['wells failing'].astype(int) - - ret_df['protocol'] = pd.Categorical(ret_df['protocol'], - categories=protocol_headings, - ordered=True) - - return ret_df - - -def run_qc(readname, savename, well, time_strs, args): - - assert len(time_strs) == 2 - filepath_before = os.path.join(args.data_directory, - f"{readname}_{time_strs[0]}") - json_file_before = f"{readname}_{time_strs[0]}" - - filepath_after = os.path.join(args.data_directory, - f"{readname}_{time_strs[1]}") - json_file_after = f"{readname}_{time_strs[1]}" - - logging.debug(f"loading {json_file_after} and {json_file_before}") - - before_trace = SyncropatchTrace(filepath_before, - json_file_before) - - after_trace = SyncropatchTrace(filepath_after, - json_file_after) - - assert before_trace.sampling_rate == after_trace.sampling_rate - - sampling_rate = before_trace.sampling_rate - - savedir = args.output_dir - if not os.path.exists(savedir): - os.makedirs(savedir) - - before_voltage = before_trace.get_voltage() - after_voltage = after_trace.get_voltage() - - # Assert that protocols are exactly the same - assert np.all(before_voltage == after_voltage) - - voltage = before_voltage - - sweeps = [0, 1] - raw_before_all = before_trace.get_trace_sweeps(sweeps) - raw_after_all = after_trace.get_trace_sweeps(sweeps) - - selected_wells = [] - - hergqc = hERGQC(sampling_rate=sampling_rate, - # plot_dir=plot_dir, - voltage=before_voltage) - - plot_dir = os.path.join(savedir, "debug", f"debug_{well}_{savename}") - - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - - qc_before = before_trace.get_onboard_QC_values() - qc_after = after_trace.get_onboard_QC_values() - - # Check if any cell first! - if (None in qc_before[well][0]) or (None in qc_after[well][0]): - no_cell = True - return False, pd.DataFrame() - - else: - no_cell = False - - nsweeps = before_trace.NofSweeps - assert after_trace.NofSweeps == nsweeps - - before_currents_corrected = np.empty((nsweeps, before_trace.NofSamples)) - after_currents_corrected = np.empty((nsweeps, after_trace.NofSamples)) - - before_currents = np.empty((nsweeps, before_trace.NofSamples)) - after_currents = np.empty((nsweeps, after_trace.NofSamples)) - - # Get ramp times from protocol description - voltage_protocol = VoltageProtocol.from_voltage_trace(voltage, - before_trace.get_times()) - - #  Find start of leak section - desc = voltage_protocol.get_all_sections() - ramp_locs = np.argwhere(desc[:, 2] != desc[:, 3]).flatten() - tstart = desc[ramp_locs[0], 0] - tend = voltage_protocol.get_ramps()[0][1] - - times = before_trace.get_times() - - ramp_bounds = [np.argmax(times > tstart), np.argmax(times > tend)] - - assert after_trace.NofSamples == before_trace.NofSamples - - for sweep in range(nsweeps): - before_raw = np.array(raw_before_all[well])[sweep, :] - after_raw = np.array(raw_after_all[well])[sweep, :] - - before_params1, before_leak = fit_linear_leak(before_raw, - voltage, - times, - *ramp_bounds, - save_fname=f"{well}-sweep{sweep}-before.png", - output_dir=savedir) - - after_params1, after_leak = fit_linear_leak(after_raw, - voltage, - times, - *ramp_bounds, - save_fname=f"{well}-sweep{sweep}-after.png", - output_dir=savedir) - - before_currents_corrected[sweep, :] = before_raw - before_leak - after_currents_corrected[sweep, :] = after_raw - after_leak - - before_currents[sweep, :] = before_raw - after_currents[sweep, :] = after_raw - - logging.info(f"{well} {savename}\n----------") - logging.info(f"sampling_rate is {sampling_rate}") - - voltage_steps = [tend - for tstart, tend, vstart, vend in - voltage_protocol.get_all_sections() if vend == vstart] - - # Run QC with raw currents - passed, QC = hergqc.run_qc(voltage_steps, times, - before_currents_corrected, - after_currents_corrected, - np.array(qc_before[well])[0, :], - np.array(qc_after[well])[0, :], nsweeps) - - QC = list(QC) - df_rows = [[well] + list(QC)] - - selected = np.all(QC) and not no_cell - if selected: - selected_wells.append(well) - - # Save subtracted current in csv file - header = "\"current\"" - - for i in range(nsweeps): - - savepath = os.path.join(savedir, - f"{args.saveID}-{savename}-{well}-sweep{i}.csv") - subtracted_current = before_currents_corrected[i, :] - after_currents_corrected[i, :] - - if args.output_traces: - if not os.path.exists(savedir): - os.makedirs(savedir) - - np.savetxt(savepath, subtracted_current, delimiter=',', - comments='', header=header) - - column_labels = ['well', 'qc1.rseal', 'qc1.cm', 'qc1.rseries', 'qc2.raw', - 'qc2.subtracted', 'qc3.raw', 'qc3.E4031', 'qc3.subtracted', - 'qc4.rseal', 'qc4.cm', 'qc4.rseries', 'qc5.staircase', - 'qc5.1.staircase', 'qc6.subtracted', 'qc6.1.subtracted', - 'qc6.2.subtracted'] - - df = pd.DataFrame(np.array(df_rows), columns=column_labels) - - missing_wells_dfs = [] - # Add onboard qc to dataframe - for well in args.wells: - if well not in df['well'].values: - onboard_qc_df = pd.DataFrame([[well] + [False for col in - list(df)[1:]]], - columns=list(df)) - missing_wells_dfs.append(onboard_qc_df) - df = pd.concat([df] + missing_wells_dfs, ignore_index=True) - - df['protocol'] = savename - - return passed, df - - -def extract_protocol(readname, savename, well, time_strs, args): - logging.info(f"extracting {savename}") - savedir = args.output_dir - saveID = args.saveID - - traces_dir = os.path.join(savedir, 'traces') - - if not os.path.exists(traces_dir): - try: - os.makedirs(traces_dir) - except FileExistsError: - pass - - row_dict = {} - - subtraction_plots_dir = os.path.join(savedir, 'subtraction_plots') - - if not os.path.isdir(subtraction_plots_dir): - try: - os.makedirs(subtraction_plots_dir) - except FileExistsError: - pass - - logging.info(f"Exporting {readname} as {savename}") - - filepath_before = os.path.join(args.data_directory, - f"{readname}_{time_strs[0]}") - filepath_after = os.path.join(args.data_directory, - f"{readname}_{time_strs[1]}") - json_file_before = f"{readname}_{time_strs[0]}" - json_file_after = f"{readname}_{time_strs[1]}" - before_trace = SyncropatchTrace(filepath_before, - json_file_before) - after_trace = SyncropatchTrace(filepath_after, - json_file_after) - - voltage_protocol = before_trace.get_voltage_protocol() - times = before_trace.get_times() - voltages = before_trace.get_voltage() - - #  Find start of leak section - desc = voltage_protocol.get_all_sections() - ramp_bounds = detect_ramp_bounds(times, desc) - tstart, tend = ramp_bounds - - nsweeps_before = before_trace.NofSweeps = 2 - nsweeps_after = after_trace.NofSweeps = 2 - - assert nsweeps_before == nsweeps_after - - # Time points - times_before = before_trace.get_times() - times_after = after_trace.get_times() - - try: - assert all(np.abs(times_before - times_after) < 1e-8) - except Exception as exc: - logging.warning(f"Exception thrown when handling {savename}: ", str(exc)) - return - - header = "\"current\"" - - qc_before = before_trace.get_onboard_QC_values() - qc_after = after_trace.get_onboard_QC_values() - qc_vals_all = before_trace.get_onboard_QC_values() - - voltage_before = before_trace.get_voltage() - voltage_after = after_trace.get_voltage() - - assert len(voltage_before) == len(voltage_after) - assert len(voltage_before) == len(times_before) - assert len(voltage_after) == len(times_after) - - voltage = voltage_before - - voltage_df = pd.DataFrame(np.vstack((times_before.flatten(), - voltage.flatten())).T, - columns=['time', 'voltage']) - - if not os.path.exists(os.path.join(traces_dir, - f"{saveID}-{savename}-voltages.csv")): - voltage_df.to_csv(os.path.join(traces_dir, - f"{saveID}-{savename}-voltages.csv")) - - np.savetxt(os.path.join(traces_dir, f"{saveID}-{savename}-times.csv"), - times_before) - - # plot subtraction - fig = plt.figure(figsize=args.figsize, layout='constrained') - - reversal_plot_dir = os.path.join(savedir, 'reversal_plots') - - rows = [] - - - before_current = before_trace.get_trace_sweeps()[well] - after_current = after_trace.get_trace_sweeps()[well] - - before_leak_currents = [] - after_leak_currents = [] - - out_dir = os.path.join(savedir, - f"{saveID}-{savename}-leak_fit-before") - - for sweep in range(before_current.shape[0]): - row_dict = { - 'well': well, - 'sweep': sweep, - 'protocol': savename - } - - qc_vals = qc_vals_all[well][sweep] - if qc_vals is None: - continue - if len(qc_vals) == 0: - continue - - row_dict['Rseal'] = qc_vals[0] - row_dict['Cm'] = qc_vals[1] - row_dict['Rseries'] = qc_vals[2] - - before_params, before_leak = fit_linear_leak(before_current[sweep, :], - voltages, times, - *ramp_bounds, - output_dir=out_dir, - save_fname=f"{well}_sweep{sweep}.png" - ) - - before_leak_currents.append(before_leak) - - out_dir = os.path.join(savedir, - f"{saveID}-{savename}-leak_fit-after") - # Convert linear regression parameters into conductance and reversal - row_dict['gleak_before'] = before_params[1] - row_dict['E_leak_before'] = -before_params[0] / before_params[1] - - after_params, after_leak = fit_linear_leak(after_current[sweep, :], - voltages, times, - *ramp_bounds, - save_fname=f"{well}_sweep{sweep}.png", - output_dir=out_dir) - - after_leak_currents.append(after_leak) - - # Convert linear regression parameters into conductance and reversal - row_dict['gleak_after'] = after_params[1] - row_dict['E_leak_after'] = -after_params[0] / after_params[1] - - subtracted_trace = before_current[sweep, :] - before_leak\ - - (after_current[sweep, :] - after_leak) - after_corrected = after_current[sweep, :] - after_leak - before_corrected = before_current[sweep, :] - before_leak - - E_rev_before = infer_reversal_potential(before_corrected, times, - desc, voltages, plot=True, - output_path=os.path.join(reversal_plot_dir, - f"{well}_{savename}_sweep{sweep}_before"), - known_Erev=args.Erev) - - E_rev_after = infer_reversal_potential(after_corrected, times, - desc, voltages, - plot=True, - output_path=os.path.join(reversal_plot_dir, - f"{well}_{savename}_sweep{sweep}_after"), - known_Erev=args.Erev) - - E_rev = infer_reversal_potential(subtracted_trace, times, desc, - voltages, plot=True, - output_path=os.path.join(reversal_plot_dir, - f"{well}_{savename}_sweep{sweep}_subtracted"), - known_Erev=args.Erev) - - row_dict['R_leftover'] =\ - np.sqrt(np.sum((after_corrected)**2)/(np.sum(before_corrected**2))) - - QC_R_leftover = np.all(row_dict['R_leftover'] < 0.5) - row_dict['QC.R_leftover'] = QC_R_leftover - - row_dict['E_rev'] = E_rev - row_dict['E_rev_before'] = E_rev_before - row_dict['E_rev_after'] = E_rev_after - - row_dict['QC.Erev'] = E_rev < -50 and E_rev > -120 - - # Check QC6 for each protocol (not just the staircase) - plot_dir = os.path.join(savedir, 'debug') - - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - - hergqc = hERGQC(sampling_rate=before_trace.sampling_rate, - plot_dir=plot_dir, - n_sweeps=before_trace.NofSweeps) - - times = before_trace.get_times() - voltage = before_trace.get_voltage() - voltage_protocol = before_trace.get_voltage_protocol() - - voltage_steps = [tstart - for tstart, tend, vstart, vend in - voltage_protocol.get_all_sections() if vend == vstart] - - current = hergqc.filter_capacitive_spikes(before_corrected - after_corrected, - times, voltage_steps) - - if args.output_traces: - out_fname = os.path.join(traces_dir, - f"{saveID}-{savename}-{well}-sweep{sweep}-subtracted.csv") - - np.savetxt(out_fname, subtracted_trace.flatten()) - rows.append(row_dict) - - extract_df = pd.DataFrame.from_dict(rows) - logging.debug(extract_df) - - times = before_trace.get_times() - voltages = before_trace.get_voltage() - - before_current_all = before_trace.get_trace_sweeps() - after_current_all = after_trace.get_trace_sweeps() - - before_current = before_current_all[well] - after_current = after_current_all[well] - - sub_df = extract_df[extract_df.well == well] - - if len(sub_df.index) == 0: - return pd.DataFrame() - - sweeps = sorted(list(sub_df.sweep.unique())) - - do_subtraction_plot(fig, times, sweeps, before_current, after_current, - voltages, ramp_bounds, well=well, - protocol=savename) - - fig.savefig(os.path.join(subtraction_plots_dir, - f"{saveID}-{savename}-{well}-sweep{sweep}-subtraction")) - fig.clf() - plt.close(fig) - - protocol_dir = os.path.join(traces_dir, 'protocols') - if not os.path.exists(protocol_dir): - try: - os.makedirs(protocol_dir) - except FileExistsError: - pass - - # extract protocol - protocol = before_trace.get_voltage_protocol() - protocol.export_txt(os.path.join(protocol_dir, - f"{saveID}-{savename}.txt")) - - json_protocol = before_trace.get_voltage_protocol_json() - - with open(os.path.join(protocol_dir, f"{saveID}-{savename}.json"), 'w') as fout: - json.dump(json_protocol, fout) - - return extract_df - - -if __name__ == '__main__': - main() From 64a90d8b74a3326d6fe5f50b28ae1381396aa556 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Wed, 27 Nov 2024 17:27:59 +0000 Subject: [PATCH 22/33] #21 Run qc1 on all test data --- pcpostprocess/hergQC.py | 2 +- pcpostprocess/scripts/run_herg_qc.py | 14 +- tests/test_herg_qc.py | 201 +++++++++++++++++++++++---- 3 files changed, 187 insertions(+), 30 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index b17fd2d..f07387d 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -1,6 +1,5 @@ import logging import os - from collections import OrderedDict import matplotlib.pyplot as plt @@ -9,6 +8,7 @@ NOISE_LEN = 200 + class QCDict: labels = [ diff --git a/pcpostprocess/scripts/run_herg_qc.py b/pcpostprocess/scripts/run_herg_qc.py index 2a400f0..75c313b 100644 --- a/pcpostprocess/scripts/run_herg_qc.py +++ b/pcpostprocess/scripts/run_herg_qc.py @@ -921,11 +921,15 @@ def run_qc_for_protocol(readname, savename, time_strs, args): voltage_protocol.get_all_sections() if vend == vstart] # Run QC with raw currents - QC = hergqc.run_qc(voltage_steps, times, - before_currents_corrected, - after_currents_corrected, - np.array(qc_before[well])[0, :], - np.array(qc_after[well])[0, :], nsweeps) + QC = hergqc.run_qc( + voltage_steps, + times, + before_currents_corrected, + after_currents_corrected, + np.array(qc_before[well])[0, :], + np.array(qc_after[well])[0, :], + nsweeps, + ) df_rows.append([well] + QC.passed_list()) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index dec3876..4ad3190 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -1,13 +1,12 @@ +import copy import logging import os -import copy -import string import unittest import numpy as np from syncropatch_export.trace import Trace -from pcpostprocess.hergQC import hERGQC, NOISE_LEN +from pcpostprocess.hergQC import NOISE_LEN, hERGQC class TestHergQC(unittest.TestCase): @@ -108,28 +107,182 @@ def passed(result): f"QC1: {rseal}, {cm}, {rseries}", ) - # Test on select data - test_matrix = [ - ("A01", True), - ("A02", True), - ("A03", True), - ("A04", True), - ("A05", True), - ("A10", False), - ("A12", False), - ("A13", False), - ("A16", False), - ("A19", False), - ] + # Test on data + test_wells_before = { + 'A01': True, 'A02': True, 'A03': True, 'A04': True, 'A05': True, + 'A06': True, 'A07': True, 'A08': True, 'A09': True, 'A10': False, + 'A11': True, 'A12': False, 'A13': False, 'A14': True, 'A15': True, + 'A16': False, 'A17': True, 'A18': True, 'A19': False, 'A20': False, + 'A21': True, 'A22': True, 'A23': True, 'A24': False, 'B01': True, + 'B02': True, 'B03': True, 'B04': True, 'B05': False, 'B06': True, + 'B07': False, 'B08': True, 'B09': True, 'B10': True, 'B11': False, + 'B12': False, 'B13': False, 'B14': True, 'B15': False, 'B16': True, + 'B17': True, 'B18': True, 'B19': False, 'B20': True, 'B21': False, + 'B22': True, 'B23': False, 'B24': True, 'C01': True, 'C02': False, + 'C03': True, 'C04': False, 'C05': True, 'C06': True, 'C07': False, + 'C08': True, 'C09': False, 'C10': True, 'C11': False, 'C12': False, + 'C13': True, 'C14': False, 'C15': True, 'C16': True, 'C17': True, + 'C18': False, 'C19': False, 'C20': False, 'C21': True, 'C22': True, + 'C23': False, 'C24': True, 'D01': True, 'D02': False, 'D03': False, + 'D04': True, 'D05': False, 'D06': True, 'D07': True, 'D08': True, + 'D09': False, 'D10': False, 'D11': True, 'D12': True, 'D13': True, + 'D14': False, 'D15': False, 'D16': False, 'D17': True, 'D18': True, + 'D19': False, 'D20': True, 'D21': False, 'D22': True, 'D23': True, + 'D24': True, 'E01': True, 'E02': True, 'E03': True, 'E04': False, + 'E05': True, 'E06': False, 'E07': False, 'E08': True, 'E09': True, + 'E10': False, 'E11': False, 'E12': True, 'E13': True, 'E14': False, + 'E15': False, 'E16': False, 'E17': False, 'E18': True, 'E19': False, + 'E20': True, 'E21': True, 'E22': False, 'E23': False, 'E24': True, + 'F01': False, 'F02': True, 'F03': False, 'F04': False, 'F05': False, + 'F06': True, 'F07': False, 'F08': True, 'F09': False, 'F10': True, + 'F11': True, 'F12': False, 'F13': False, 'F14': False, 'F15': False, + 'F16': True, 'F17': True, 'F18': False, 'F19': False, 'F20': False, + 'F21': False, 'F22': True, 'F23': True, 'F24': False, 'G01': True, + 'G02': True, 'G03': True, 'G04': True, 'G05': True, 'G06': False, + 'G07': True, 'G08': True, 'G09': False, 'G10': True, 'G11': True, + 'G12': False, 'G13': False, 'G14': False, 'G15': True, 'G16': False, + 'G17': False, 'G18': True, 'G19': True, 'G20': False, 'G21': False, + 'G22': True, 'G23': False, 'G24': False, 'H01': False, 'H02': False, + 'H03': False, 'H04': True, 'H05': True, 'H06': False, 'H07': False, + 'H08': False, 'H09': True, 'H10': False, 'H11': False, 'H12': True, + 'H13': False, 'H14': False, 'H15': False, 'H16': False, 'H17': True, + 'H18': True, 'H19': False, 'H20': True, 'H21': False, 'H22': True, + 'H23': False, 'H24': False, 'I01': False, 'I02': True, 'I03': True, + 'I04': False, 'I05': False, 'I06': False, 'I07': False, 'I08': False, + 'I09': True, 'I10': False, 'I11': False, 'I12': False, 'I13': True, + 'I14': True, 'I15': True, 'I16': False, 'I17': False, 'I18': True, + 'I19': True, 'I20': True, 'I21': False, 'I22': True, 'I23': True, + 'I24': True, 'J01': True, 'J02': True, 'J03': True, 'J04': True, + 'J05': True, 'J06': True, 'J07': False, 'J08': True, 'J09': True, + 'J10': False, 'J11': True, 'J12': True, 'J13': True, 'J14': True, + 'J15': True, 'J16': False, 'J17': False, 'J18': True, 'J19': False, + 'J20': True, 'J21': False, 'J22': True, 'J23': True, 'J24': False, + 'K01': True, 'K02': False, 'K03': False, 'K04': True, 'K05': True, + 'K06': False, 'K07': False, 'K08': True, 'K09': True, 'K10': True, + 'K11': False, 'K12': False, 'K13': True, 'K14': True, 'K15': True, + 'K16': False, 'K17': False, 'K18': True, 'K19': True, 'K20': False, + 'K21': True, 'K22': False, 'K23': True, 'K24': False, 'L01': False, + 'L02': False, 'L03': True, 'L04': False, 'L05': False, 'L06': True, + 'L07': True, 'L08': False, 'L09': True, 'L10': False, 'L11': False, + 'L12': True, 'L13': False, 'L14': True, 'L15': True, 'L16': False, + 'L17': False, 'L18': False, 'L19': True, 'L20': True, 'L21': True, + 'L22': True, 'L23': True, 'L24': False, 'M01': False, 'M02': True, + 'M03': True, 'M04': False, 'M05': True, 'M06': False, 'M07': True, + 'M08': True, 'M09': False, 'M10': True, 'M11': True, 'M12': False, + 'M13': True, 'M14': False, 'M15': False, 'M16': False, 'M17': True, + 'M18': True, 'M19': False, 'M20': False, 'M21': False, 'M22': True, + 'M23': True, 'M24': True, 'N01': True, 'N02': True, 'N03': False, + 'N04': False, 'N05': True, 'N06': False, 'N07': True, 'N08': False, + 'N09': True, 'N10': True, 'N11': False, 'N12': True, 'N13': False, + 'N14': False, 'N15': True, 'N16': False, 'N17': True, 'N18': False, + 'N19': True, 'N20': True, 'N21': False, 'N22': True, 'N23': True, + 'N24': False, 'O01': False, 'O02': False, 'O03': False, 'O04': True, + 'O05': False, 'O06': True, 'O07': False, 'O08': True, 'O09': True, + 'O10': False, 'O11': False, 'O12': True, 'O13': True, 'O14': True, + 'O15': True, 'O16': True, 'O17': False, 'O18': True, 'O19': False, + 'O20': True, 'O21': True, 'O22': False, 'O23': True, 'O24': False, + 'P01': False, 'P02': True, 'P03': False, 'P04': True, 'P05': True, + 'P06': False, 'P07': False, 'P08': False, 'P09': False, 'P10': True, + 'P11': True, 'P12': False, 'P13': False, 'P14': False, 'P15': False, + 'P16': False, 'P17': False, 'P18': False, 'P19': True, 'P20': True, + 'P21': False, 'P22': False, 'P23': True, 'P24': False + } + + for well in test_wells_before: + qc_vals_before = np.array(self.qc_vals_before[well])[0, :] + self.assertEqual( + passed(hergqc.qc1(*qc_vals_before)), + test_wells_before[well], + f"QC1: {well} (before) {qc_vals_before}", + ) - for well, expected in test_matrix: - with self.subTest(well): - qc_vals_before = np.array(self.qc_vals_before[well])[0, :] - self.assertEqual( - passed(hergqc.qc1(*qc_vals_before)), - expected, - f"QC1: {well} (before) {qc_vals_before}", - ) + test_wells_after = { + 'A01': True, 'A02': True, 'A03': True, 'A04': True, 'A05': True, + 'A06': False, 'A07': True, 'A08': False, 'A09': True, 'A10': False, + 'A11': True, 'A12': False, 'A13': False, 'A14': True, 'A15': True, + 'A16': False, 'A17': True, 'A18': True, 'A19': False, 'A20': False, + 'A21': True, 'A22': True, 'A23': True, 'A24': False, 'B01': True, + 'B02': False, 'B03': True, 'B04': True, 'B05': False, 'B06': True, + 'B07': False, 'B08': True, 'B09': True, 'B10': True, 'B11': False, + 'B12': False, 'B13': False, 'B14': True, 'B15': False, 'B16': True, + 'B17': True, 'B18': True, 'B19': False, 'B20': True, 'B21': False, + 'B22': True, 'B23': False, 'B24': True, 'C01': False, 'C02': False, + 'C03': True, 'C04': False, 'C05': True, 'C06': True, 'C07': False, + 'C08': True, 'C09': False, 'C10': True, 'C11': False, 'C12': False, + 'C13': True, 'C14': False, 'C15': True, 'C16': True, 'C17': True, + 'C18': False, 'C19': True, 'C20': False, 'C21': True, 'C22': False, + 'C23': False, 'C24': True, 'D01': True, 'D02': True, 'D03': False, + 'D04': True, 'D05': False, 'D06': True, 'D07': True, 'D08': True, + 'D09': False, 'D10': False, 'D11': True, 'D12': True, 'D13': True, + 'D14': False, 'D15': False, 'D16': True, 'D17': True, 'D18': True, + 'D19': False, 'D20': True, 'D21': False, 'D22': True, 'D23': True, + 'D24': True, 'E01': False, 'E02': True, 'E03': False, 'E04': False, + 'E05': True, 'E06': False, 'E07': False, 'E08': True, 'E09': False, + 'E10': False, 'E11': False, 'E12': True, 'E13': True, 'E14': False, + 'E15': False, 'E16': False, 'E17': False, 'E18': True, 'E19': False, + 'E20': False, 'E21': True, 'E22': False, 'E23': False, 'E24': False, + 'F01': False, 'F02': True, 'F03': False, 'F04': False, 'F05': True, + 'F06': True, 'F07': False, 'F08': True, 'F09': False, 'F10': True, + 'F11': True, 'F12': False, 'F13': False, 'F14': False, 'F15': False, + 'F16': False, 'F17': True, 'F18': False, 'F19': False, 'F20': False, + 'F21': False, 'F22': True, 'F23': True, 'F24': False, 'G01': True, + 'G02': True, 'G03': True, 'G04': True, 'G05': True, 'G06': False, + 'G07': True, 'G08': False, 'G09': False, 'G10': True, 'G11': True, + 'G12': False, 'G13': False, 'G14': False, 'G15': False, 'G16': False, + 'G17': False, 'G18': True, 'G19': True, 'G20': False, 'G21': False, + 'G22': True, 'G23': False, 'G24': False, 'H01': False, 'H02': False, + 'H03': False, 'H04': False, 'H05': True, 'H06': False, 'H07': False, + 'H08': False, 'H09': False, 'H10': False, 'H11': False, 'H12': True, + 'H13': False, 'H14': False, 'H15': False, 'H16': False, 'H17': False, + 'H18': False, 'H19': False, 'H20': False, 'H21': False, 'H22': True, + 'H23': False, 'H24': False, 'I01': False, 'I02': True, 'I03': True, + 'I04': False, 'I05': False, 'I06': False, 'I07': False, 'I08': False, + 'I09': True, 'I10': False, 'I11': True, 'I12': False, 'I13': True, + 'I14': True, 'I15': True, 'I16': False, 'I17': False, 'I18': True, + 'I19': True, 'I20': False, 'I21': False, 'I22': True, 'I23': True, + 'I24': True, 'J01': True, 'J02': True, 'J03': False, 'J04': True, + 'J05': True, 'J06': True, 'J07': False, 'J08': True, 'J09': True, + 'J10': False, 'J11': True, 'J12': True, 'J13': True, 'J14': False, + 'J15': True, 'J16': False, 'J17': False, 'J18': True, 'J19': False, + 'J20': True, 'J21': False, 'J22': True, 'J23': False, 'J24': False, + 'K01': True, 'K02': False, 'K03': False, 'K04': True, 'K05': True, + 'K06': False, 'K07': False, 'K08': True, 'K09': True, 'K10': True, + 'K11': True, 'K12': True, 'K13': True, 'K14': False, 'K15': True, + 'K16': False, 'K17': True, 'K18': True, 'K19': True, 'K20': False, + 'K21': True, 'K22': False, 'K23': False, 'K24': False, 'L01': False, + 'L02': False, 'L03': False, 'L04': False, 'L05': False, 'L06': True, + 'L07': True, 'L08': False, 'L09': True, 'L10': False, 'L11': False, + 'L12': False, 'L13': False, 'L14': True, 'L15': True, 'L16': False, + 'L17': False, 'L18': False, 'L19': True, 'L20': True, 'L21': True, + 'L22': True, 'L23': True, 'L24': False, 'M01': False, 'M02': False, + 'M03': True, 'M04': False, 'M05': False, 'M06': False, 'M07': True, + 'M08': False, 'M09': False, 'M10': True, 'M11': True, 'M12': False, + 'M13': False, 'M14': False, 'M15': False, 'M16': True, 'M17': True, + 'M18': True, 'M19': False, 'M20': False, 'M21': False, 'M22': True, + 'M23': True, 'M24': True, 'N01': True, 'N02': True, 'N03': True, + 'N04': False, 'N05': True, 'N06': False, 'N07': True, 'N08': False, + 'N09': True, 'N10': True, 'N11': False, 'N12': True, 'N13': False, + 'N14': False, 'N15': True, 'N16': False, 'N17': True, 'N18': False, + 'N19': False, 'N20': True, 'N21': False, 'N22': True, 'N23': True, + 'N24': False, 'O01': False, 'O02': False, 'O03': False, 'O04': True, + 'O05': False, 'O06': True, 'O07': False, 'O08': False, 'O09': True, + 'O10': False, 'O11': False, 'O12': True, 'O13': True, 'O14': False, + 'O15': False, 'O16': True, 'O17': False, 'O18': True, 'O19': False, + 'O20': True, 'O21': True, 'O22': False, 'O23': True, 'O24': False, + 'P01': False, 'P02': True, 'P03': False, 'P04': True, 'P05': True, + 'P06': False, 'P07': False, 'P08': False, 'P09': False, 'P10': True, + 'P11': False, 'P12': False, 'P13': False, 'P14': False, 'P15': False, + 'P16': False, 'P17': False, 'P18': False, 'P19': False, 'P20': True, + 'P21': False, 'P22': False, 'P23': True, 'P24': True + } + + for well in test_wells_after: + qc_vals_after = np.array(self.qc_vals_after[well])[0, :] + self.assertEqual( + passed(hergqc.qc1(*qc_vals_after)), + test_wells_after[well], + f"QC1: {well} (after) {qc_vals_after}", + ) def test_qc2(self): hergqc = copy.deepcopy(self.hergqc) From cf9aa1ca153c7e5d996e7f58eb95b36aa59eaa61 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Wed, 27 Nov 2024 18:38:36 +0000 Subject: [PATCH 23/33] #21 Run qc2 on all test data --- tests/test_herg_qc.py | 318 ++++++++++++++++-------------------------- 1 file changed, 124 insertions(+), 194 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 4ad3190..c9cf864 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -1,6 +1,7 @@ import copy import logging import os +import string import unittest import numpy as np @@ -9,6 +10,10 @@ from pcpostprocess.hergQC import NOISE_LEN, hERGQC +def all_passed(result): + return all([x for x, _ in result]) + + class TestHergQC(unittest.TestCase): def setUp(self): @@ -32,6 +37,13 @@ def setUp(self): self.voltage = trace_before.get_voltage() self.times = trace_after.get_times() + self.n_sweeps = 2 + + self.all_wells = [ + row + str(i).zfill(2) + for row in string.ascii_uppercase[:16] + for i in range(1, 25) + ] sampling_rate = int(1.0 / (self.times[1] - self.times[0])) # in kHz @@ -58,8 +70,6 @@ def test_qc_inputs(self): self.assertTrue(np.all(np.isfinite(self.times))) def test_qc1(self): - def passed(result): - return all([x for x, _ in result]) hergqc = copy.deepcopy(self.hergqc) plot_dir = os.path.join(hergqc.plot_dir, "test_qc1") @@ -102,185 +112,73 @@ def passed(result): for (rseal, cm, rseries), expected in test_matrix: self.assertEqual( - passed(hergqc.qc1(rseal, cm, rseries)), + all_passed(hergqc.qc1(rseal, cm, rseries)), expected, f"QC1: {rseal}, {cm}, {rseries}", ) # Test on data - test_wells_before = { - 'A01': True, 'A02': True, 'A03': True, 'A04': True, 'A05': True, - 'A06': True, 'A07': True, 'A08': True, 'A09': True, 'A10': False, - 'A11': True, 'A12': False, 'A13': False, 'A14': True, 'A15': True, - 'A16': False, 'A17': True, 'A18': True, 'A19': False, 'A20': False, - 'A21': True, 'A22': True, 'A23': True, 'A24': False, 'B01': True, - 'B02': True, 'B03': True, 'B04': True, 'B05': False, 'B06': True, - 'B07': False, 'B08': True, 'B09': True, 'B10': True, 'B11': False, - 'B12': False, 'B13': False, 'B14': True, 'B15': False, 'B16': True, - 'B17': True, 'B18': True, 'B19': False, 'B20': True, 'B21': False, - 'B22': True, 'B23': False, 'B24': True, 'C01': True, 'C02': False, - 'C03': True, 'C04': False, 'C05': True, 'C06': True, 'C07': False, - 'C08': True, 'C09': False, 'C10': True, 'C11': False, 'C12': False, - 'C13': True, 'C14': False, 'C15': True, 'C16': True, 'C17': True, - 'C18': False, 'C19': False, 'C20': False, 'C21': True, 'C22': True, - 'C23': False, 'C24': True, 'D01': True, 'D02': False, 'D03': False, - 'D04': True, 'D05': False, 'D06': True, 'D07': True, 'D08': True, - 'D09': False, 'D10': False, 'D11': True, 'D12': True, 'D13': True, - 'D14': False, 'D15': False, 'D16': False, 'D17': True, 'D18': True, - 'D19': False, 'D20': True, 'D21': False, 'D22': True, 'D23': True, - 'D24': True, 'E01': True, 'E02': True, 'E03': True, 'E04': False, - 'E05': True, 'E06': False, 'E07': False, 'E08': True, 'E09': True, - 'E10': False, 'E11': False, 'E12': True, 'E13': True, 'E14': False, - 'E15': False, 'E16': False, 'E17': False, 'E18': True, 'E19': False, - 'E20': True, 'E21': True, 'E22': False, 'E23': False, 'E24': True, - 'F01': False, 'F02': True, 'F03': False, 'F04': False, 'F05': False, - 'F06': True, 'F07': False, 'F08': True, 'F09': False, 'F10': True, - 'F11': True, 'F12': False, 'F13': False, 'F14': False, 'F15': False, - 'F16': True, 'F17': True, 'F18': False, 'F19': False, 'F20': False, - 'F21': False, 'F22': True, 'F23': True, 'F24': False, 'G01': True, - 'G02': True, 'G03': True, 'G04': True, 'G05': True, 'G06': False, - 'G07': True, 'G08': True, 'G09': False, 'G10': True, 'G11': True, - 'G12': False, 'G13': False, 'G14': False, 'G15': True, 'G16': False, - 'G17': False, 'G18': True, 'G19': True, 'G20': False, 'G21': False, - 'G22': True, 'G23': False, 'G24': False, 'H01': False, 'H02': False, - 'H03': False, 'H04': True, 'H05': True, 'H06': False, 'H07': False, - 'H08': False, 'H09': True, 'H10': False, 'H11': False, 'H12': True, - 'H13': False, 'H14': False, 'H15': False, 'H16': False, 'H17': True, - 'H18': True, 'H19': False, 'H20': True, 'H21': False, 'H22': True, - 'H23': False, 'H24': False, 'I01': False, 'I02': True, 'I03': True, - 'I04': False, 'I05': False, 'I06': False, 'I07': False, 'I08': False, - 'I09': True, 'I10': False, 'I11': False, 'I12': False, 'I13': True, - 'I14': True, 'I15': True, 'I16': False, 'I17': False, 'I18': True, - 'I19': True, 'I20': True, 'I21': False, 'I22': True, 'I23': True, - 'I24': True, 'J01': True, 'J02': True, 'J03': True, 'J04': True, - 'J05': True, 'J06': True, 'J07': False, 'J08': True, 'J09': True, - 'J10': False, 'J11': True, 'J12': True, 'J13': True, 'J14': True, - 'J15': True, 'J16': False, 'J17': False, 'J18': True, 'J19': False, - 'J20': True, 'J21': False, 'J22': True, 'J23': True, 'J24': False, - 'K01': True, 'K02': False, 'K03': False, 'K04': True, 'K05': True, - 'K06': False, 'K07': False, 'K08': True, 'K09': True, 'K10': True, - 'K11': False, 'K12': False, 'K13': True, 'K14': True, 'K15': True, - 'K16': False, 'K17': False, 'K18': True, 'K19': True, 'K20': False, - 'K21': True, 'K22': False, 'K23': True, 'K24': False, 'L01': False, - 'L02': False, 'L03': True, 'L04': False, 'L05': False, 'L06': True, - 'L07': True, 'L08': False, 'L09': True, 'L10': False, 'L11': False, - 'L12': True, 'L13': False, 'L14': True, 'L15': True, 'L16': False, - 'L17': False, 'L18': False, 'L19': True, 'L20': True, 'L21': True, - 'L22': True, 'L23': True, 'L24': False, 'M01': False, 'M02': True, - 'M03': True, 'M04': False, 'M05': True, 'M06': False, 'M07': True, - 'M08': True, 'M09': False, 'M10': True, 'M11': True, 'M12': False, - 'M13': True, 'M14': False, 'M15': False, 'M16': False, 'M17': True, - 'M18': True, 'M19': False, 'M20': False, 'M21': False, 'M22': True, - 'M23': True, 'M24': True, 'N01': True, 'N02': True, 'N03': False, - 'N04': False, 'N05': True, 'N06': False, 'N07': True, 'N08': False, - 'N09': True, 'N10': True, 'N11': False, 'N12': True, 'N13': False, - 'N14': False, 'N15': True, 'N16': False, 'N17': True, 'N18': False, - 'N19': True, 'N20': True, 'N21': False, 'N22': True, 'N23': True, - 'N24': False, 'O01': False, 'O02': False, 'O03': False, 'O04': True, - 'O05': False, 'O06': True, 'O07': False, 'O08': True, 'O09': True, - 'O10': False, 'O11': False, 'O12': True, 'O13': True, 'O14': True, - 'O15': True, 'O16': True, 'O17': False, 'O18': True, 'O19': False, - 'O20': True, 'O21': True, 'O22': False, 'O23': True, 'O24': False, - 'P01': False, 'P02': True, 'P03': False, 'P04': True, 'P05': True, - 'P06': False, 'P07': False, 'P08': False, 'P09': False, 'P10': True, - 'P11': True, 'P12': False, 'P13': False, 'P14': False, 'P15': False, - 'P16': False, 'P17': False, 'P18': False, 'P19': True, 'P20': True, - 'P21': False, 'P22': False, 'P23': True, 'P24': False - } - - for well in test_wells_before: + failed_wells_before = [ + 'A10', 'A12', 'A13', 'A16', 'A19', 'A20', 'A24', 'B05', 'B07', 'B11', + 'B12', 'B13', 'B15', 'B19', 'B21', 'B23', 'C02', 'C04', 'C07', 'C09', + 'C11', 'C12', 'C14', 'C18', 'C19', 'C20', 'C23', 'D02', 'D03', 'D05', + 'D09', 'D10', 'D14', 'D15', 'D16', 'D19', 'D21', 'E04', 'E06', 'E07', + 'E10', 'E11', 'E14', 'E15', 'E16', 'E17', 'E19', 'E22', 'E23', 'F01', + 'F03', 'F04', 'F05', 'F07', 'F09', 'F12', 'F13', 'F14', 'F15', 'F18', + 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', 'G13', 'G14', 'G16', + 'G17', 'G20', 'G21', 'G23', 'G24', 'H01', 'H02', 'H03', 'H06', 'H07', + 'H08', 'H10', 'H11', 'H13', 'H14', 'H15', 'H16', 'H19', 'H21', 'H23', + 'H24', 'I01', 'I04', 'I05', 'I06', 'I07', 'I08', 'I10', 'I11', 'I12', + 'I16', 'I17', 'I21', 'J07', 'J10', 'J16', 'J17', 'J19', 'J21', 'J24', + 'K02', 'K03', 'K06', 'K07', 'K11', 'K12', 'K16', 'K17', 'K20', 'K22', + 'K24', 'L01', 'L02', 'L04', 'L05', 'L08', 'L10', 'L11', 'L13', 'L16', + 'L17', 'L18', 'L24', 'M01', 'M04', 'M06', 'M09', 'M12', 'M14', 'M15', + 'M16', 'M19', 'M20', 'M21', 'N03', 'N04', 'N06', 'N08', 'N11', 'N13', + 'N14', 'N16', 'N18', 'N21', 'N24', 'O01', 'O02', 'O03', 'O05', 'O07', + 'O10', 'O11', 'O17', 'O19', 'O22', 'O24', 'P01', 'P03', 'P06', 'P07', + 'P08', 'P09', 'P12', 'P13', 'P14', 'P15', 'P16', 'P17', 'P18', 'P21', + 'P22', 'P24' + ] + + for well in self.all_wells: qc_vals_before = np.array(self.qc_vals_before[well])[0, :] + ex_pass_before = well not in failed_wells_before self.assertEqual( - passed(hergqc.qc1(*qc_vals_before)), - test_wells_before[well], + all_passed(hergqc.qc1(*qc_vals_before)), + ex_pass_before, f"QC1: {well} (before) {qc_vals_before}", ) - test_wells_after = { - 'A01': True, 'A02': True, 'A03': True, 'A04': True, 'A05': True, - 'A06': False, 'A07': True, 'A08': False, 'A09': True, 'A10': False, - 'A11': True, 'A12': False, 'A13': False, 'A14': True, 'A15': True, - 'A16': False, 'A17': True, 'A18': True, 'A19': False, 'A20': False, - 'A21': True, 'A22': True, 'A23': True, 'A24': False, 'B01': True, - 'B02': False, 'B03': True, 'B04': True, 'B05': False, 'B06': True, - 'B07': False, 'B08': True, 'B09': True, 'B10': True, 'B11': False, - 'B12': False, 'B13': False, 'B14': True, 'B15': False, 'B16': True, - 'B17': True, 'B18': True, 'B19': False, 'B20': True, 'B21': False, - 'B22': True, 'B23': False, 'B24': True, 'C01': False, 'C02': False, - 'C03': True, 'C04': False, 'C05': True, 'C06': True, 'C07': False, - 'C08': True, 'C09': False, 'C10': True, 'C11': False, 'C12': False, - 'C13': True, 'C14': False, 'C15': True, 'C16': True, 'C17': True, - 'C18': False, 'C19': True, 'C20': False, 'C21': True, 'C22': False, - 'C23': False, 'C24': True, 'D01': True, 'D02': True, 'D03': False, - 'D04': True, 'D05': False, 'D06': True, 'D07': True, 'D08': True, - 'D09': False, 'D10': False, 'D11': True, 'D12': True, 'D13': True, - 'D14': False, 'D15': False, 'D16': True, 'D17': True, 'D18': True, - 'D19': False, 'D20': True, 'D21': False, 'D22': True, 'D23': True, - 'D24': True, 'E01': False, 'E02': True, 'E03': False, 'E04': False, - 'E05': True, 'E06': False, 'E07': False, 'E08': True, 'E09': False, - 'E10': False, 'E11': False, 'E12': True, 'E13': True, 'E14': False, - 'E15': False, 'E16': False, 'E17': False, 'E18': True, 'E19': False, - 'E20': False, 'E21': True, 'E22': False, 'E23': False, 'E24': False, - 'F01': False, 'F02': True, 'F03': False, 'F04': False, 'F05': True, - 'F06': True, 'F07': False, 'F08': True, 'F09': False, 'F10': True, - 'F11': True, 'F12': False, 'F13': False, 'F14': False, 'F15': False, - 'F16': False, 'F17': True, 'F18': False, 'F19': False, 'F20': False, - 'F21': False, 'F22': True, 'F23': True, 'F24': False, 'G01': True, - 'G02': True, 'G03': True, 'G04': True, 'G05': True, 'G06': False, - 'G07': True, 'G08': False, 'G09': False, 'G10': True, 'G11': True, - 'G12': False, 'G13': False, 'G14': False, 'G15': False, 'G16': False, - 'G17': False, 'G18': True, 'G19': True, 'G20': False, 'G21': False, - 'G22': True, 'G23': False, 'G24': False, 'H01': False, 'H02': False, - 'H03': False, 'H04': False, 'H05': True, 'H06': False, 'H07': False, - 'H08': False, 'H09': False, 'H10': False, 'H11': False, 'H12': True, - 'H13': False, 'H14': False, 'H15': False, 'H16': False, 'H17': False, - 'H18': False, 'H19': False, 'H20': False, 'H21': False, 'H22': True, - 'H23': False, 'H24': False, 'I01': False, 'I02': True, 'I03': True, - 'I04': False, 'I05': False, 'I06': False, 'I07': False, 'I08': False, - 'I09': True, 'I10': False, 'I11': True, 'I12': False, 'I13': True, - 'I14': True, 'I15': True, 'I16': False, 'I17': False, 'I18': True, - 'I19': True, 'I20': False, 'I21': False, 'I22': True, 'I23': True, - 'I24': True, 'J01': True, 'J02': True, 'J03': False, 'J04': True, - 'J05': True, 'J06': True, 'J07': False, 'J08': True, 'J09': True, - 'J10': False, 'J11': True, 'J12': True, 'J13': True, 'J14': False, - 'J15': True, 'J16': False, 'J17': False, 'J18': True, 'J19': False, - 'J20': True, 'J21': False, 'J22': True, 'J23': False, 'J24': False, - 'K01': True, 'K02': False, 'K03': False, 'K04': True, 'K05': True, - 'K06': False, 'K07': False, 'K08': True, 'K09': True, 'K10': True, - 'K11': True, 'K12': True, 'K13': True, 'K14': False, 'K15': True, - 'K16': False, 'K17': True, 'K18': True, 'K19': True, 'K20': False, - 'K21': True, 'K22': False, 'K23': False, 'K24': False, 'L01': False, - 'L02': False, 'L03': False, 'L04': False, 'L05': False, 'L06': True, - 'L07': True, 'L08': False, 'L09': True, 'L10': False, 'L11': False, - 'L12': False, 'L13': False, 'L14': True, 'L15': True, 'L16': False, - 'L17': False, 'L18': False, 'L19': True, 'L20': True, 'L21': True, - 'L22': True, 'L23': True, 'L24': False, 'M01': False, 'M02': False, - 'M03': True, 'M04': False, 'M05': False, 'M06': False, 'M07': True, - 'M08': False, 'M09': False, 'M10': True, 'M11': True, 'M12': False, - 'M13': False, 'M14': False, 'M15': False, 'M16': True, 'M17': True, - 'M18': True, 'M19': False, 'M20': False, 'M21': False, 'M22': True, - 'M23': True, 'M24': True, 'N01': True, 'N02': True, 'N03': True, - 'N04': False, 'N05': True, 'N06': False, 'N07': True, 'N08': False, - 'N09': True, 'N10': True, 'N11': False, 'N12': True, 'N13': False, - 'N14': False, 'N15': True, 'N16': False, 'N17': True, 'N18': False, - 'N19': False, 'N20': True, 'N21': False, 'N22': True, 'N23': True, - 'N24': False, 'O01': False, 'O02': False, 'O03': False, 'O04': True, - 'O05': False, 'O06': True, 'O07': False, 'O08': False, 'O09': True, - 'O10': False, 'O11': False, 'O12': True, 'O13': True, 'O14': False, - 'O15': False, 'O16': True, 'O17': False, 'O18': True, 'O19': False, - 'O20': True, 'O21': True, 'O22': False, 'O23': True, 'O24': False, - 'P01': False, 'P02': True, 'P03': False, 'P04': True, 'P05': True, - 'P06': False, 'P07': False, 'P08': False, 'P09': False, 'P10': True, - 'P11': False, 'P12': False, 'P13': False, 'P14': False, 'P15': False, - 'P16': False, 'P17': False, 'P18': False, 'P19': False, 'P20': True, - 'P21': False, 'P22': False, 'P23': True, 'P24': True - } - - for well in test_wells_after: + failed_wells_after = [ + 'A06', 'A08', 'A10', 'A12', 'A13', 'A16', 'A19', 'A20', 'A24', 'B02', + 'B05', 'B07', 'B11', 'B12', 'B13', 'B15', 'B19', 'B21', 'B23', 'C01', + 'C02', 'C04', 'C07', 'C09', 'C11', 'C12', 'C14', 'C18', 'C20', 'C22', + 'C23', 'D03', 'D05', 'D09', 'D10', 'D14', 'D15', 'D19', 'D21', 'E01', + 'E03', 'E04', 'E06', 'E07', 'E09', 'E10', 'E11', 'E14', 'E15', 'E16', + 'E17', 'E19', 'E20', 'E22', 'E23', 'E24', 'F01', 'F03', 'F04', 'F07', + 'F09', 'F12', 'F13', 'F14', 'F15', 'F16', 'F18', 'F19', 'F20', 'F21', + 'F24', 'G06', 'G08', 'G09', 'G12', 'G13', 'G14', 'G15', 'G16', 'G17', + 'G20', 'G21', 'G23', 'G24', 'H01', 'H02', 'H03', 'H04', 'H06', 'H07', + 'H08', 'H09', 'H10', 'H11', 'H13', 'H14', 'H15', 'H16', 'H17', 'H18', + 'H19', 'H20', 'H21', 'H23', 'H24', 'I01', 'I04', 'I05', 'I06', 'I07', + 'I08', 'I10', 'I12', 'I16', 'I17', 'I20', 'I21', 'J03', 'J07', 'J10', + 'J14', 'J16', 'J17', 'J19', 'J21', 'J23', 'J24', 'K02', 'K03', 'K06', + 'K07', 'K14', 'K16', 'K20', 'K22', 'K23', 'K24', 'L01', 'L02', 'L03', + 'L04', 'L05', 'L08', 'L10', 'L11', 'L12', 'L13', 'L16', 'L17', 'L18', + 'L24', 'M01', 'M02', 'M04', 'M05', 'M06', 'M08', 'M09', 'M12', 'M13', + 'M14', 'M15', 'M19', 'M20', 'M21', 'N04', 'N06', 'N08', 'N11', 'N13', + 'N14', 'N16', 'N18', 'N19', 'N21', 'N24', 'O01', 'O02', 'O03', 'O05', + 'O07', 'O08', 'O10', 'O11', 'O14', 'O15', 'O17', 'O19', 'O22', 'O24', + 'P01', 'P03', 'P06', 'P07', 'P08', 'P09', 'P11', 'P12', 'P13', 'P14', + 'P15', 'P16', 'P17', 'P18', 'P19', 'P21', 'P22' + ] + + for well in self.all_wells: qc_vals_after = np.array(self.qc_vals_after[well])[0, :] + ex_pass_after = well not in failed_wells_after self.assertEqual( - passed(hergqc.qc1(*qc_vals_after)), - test_wells_after[well], + all_passed(hergqc.qc1(*qc_vals_after)), + ex_pass_after, f"QC1: {well} (after) {qc_vals_after}", ) @@ -293,20 +191,50 @@ def test_qc2(self): # qc2 checks that raw and subtracted SNR are above a minimum threshold test_matrix = [ - (10, True), # snr = 8082 - (1, True), # snr = 74 - (0.601, True), # snr = 25.07 - (0.6, False), # snr = 24.98 - (0.5, False), # snr = 17 - (0.1, False), # snr = 0.5 + (10, 8082.1, True), + (1, 74.0, True), + (0.601, 25.1, True), + (0.6, 25.0, False), + (0.5, 16.8, False), + (0.1, 0.5, False), ] - for i, expected in test_matrix: + for i, ex_snr, ex_pass in test_matrix: recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [i] * 500) - result = hergqc.qc2(recording) - self.assertEqual(result[0], expected, f"({i}: {result[1]})") + pass_, snr = hergqc.qc2(recording) + self.assertAlmostEqual( + snr, ex_snr, 1, f"QC2: ({i}) {snr} != {ex_snr}") + self.assertEqual(pass_, ex_pass, f"QC2: ({i}) {pass_} != {ex_pass}") - # TODO: Test on select data + # Test on data + failed_wells_raw = ["P16"] + failed_wells_subtracted = [ + "B09", "C11", "H19", "H24", "K22", "O16", "P16" + ] + + for well in self.all_wells: + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) + + raw = [] + subtracted = [] + for i in range(self.n_sweeps): + raw.append(hergqc.qc2(before[i])) + subtracted.append(hergqc.qc2(before[i] - after[i])) + + ex_pass_raw = well not in failed_wells_raw + self.assertEqual( + all_passed(raw), + ex_pass_raw, + f"QC2: {well} (raw) {raw}", + ) + + ex_pass_subtracted = well not in failed_wells_subtracted + self.assertEqual( + all_passed(subtracted), + ex_pass_subtracted, + f"QC2: {well} (subtracted) {subtracted}", + ) def test_qc3(self): hergqc = copy.deepcopy(self.hergqc) @@ -331,7 +259,8 @@ def test_qc3(self): recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) for i, expected in test_matrix: - recording2 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40 + i] * 500) + recording2 = np.asarray( + [0, 0.1] * (NOISE_LEN // 2) + [40 + i] * 500) result = hergqc.qc3(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -345,16 +274,14 @@ def test_qc3(self): recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) for i, expected in test_matrix: - recording2 = np.asarray([0, 0.1 * i] * (NOISE_LEN // 2) + [40] * 500) + recording2 = np.asarray( + [0, 0.1 * i] * (NOISE_LEN // 2) + [40] * 500) result = hergqc.qc3(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") # TODO: Test on select data def test_qc4(self): - def passed(result): - return all([x for x, _ in result]) - hergqc = copy.deepcopy(self.hergqc) plot_dir = os.path.join(hergqc.plot_dir, "test_qc4") if not os.path.exists(plot_dir): @@ -379,14 +306,14 @@ def passed(result): for i, expected in test_matrix: rseals = [r_lo, i * r_lo] self.assertEqual( - passed(hergqc.qc4(rseals, cms, rseriess)), + all_passed(hergqc.qc4(rseals, cms, rseriess)), expected, f"({i}: {rseals}, {cms}, {rseriess})", ) rseals = [r_hi, i * r_hi] self.assertEqual( - passed(hergqc.qc4(rseals, cms, rseriess)), + all_passed(hergqc.qc4(rseals, cms, rseriess)), expected, f"({i}: {rseals}, {cms}, {rseriess})", ) @@ -405,14 +332,14 @@ def passed(result): for i, expected in test_matrix: cms = [c_lo, i * c_lo] self.assertEqual( - passed(hergqc.qc4(rseals, cms, rseriess)), + all_passed(hergqc.qc4(rseals, cms, rseriess)), expected, f"({i}: {rseals}, {cms}, {rseriess})", ) cms = [c_hi, i * c_hi] self.assertEqual( - passed(hergqc.qc4(rseals, cms, rseriess)), + all_passed(hergqc.qc4(rseals, cms, rseriess)), expected, f"({i}: {rseals}, {cms}, {rseriess})", ) @@ -431,14 +358,14 @@ def passed(result): for i, expected in test_matrix: rseriess = [r_lo, i * r_lo] self.assertEqual( - passed(hergqc.qc4(rseals, cms, rseriess)), + all_passed(hergqc.qc4(rseals, cms, rseriess)), expected, f"({i}: {rseals}, {cms}, {rseriess})", ) rseriess = [r_hi, i * r_hi] self.assertEqual( - passed(hergqc.qc4(rseals, cms, rseriess)), + all_passed(hergqc.qc4(rseals, cms, rseriess)), expected, f"({i}: {rseals}, {cms}, {rseriess})", ) @@ -466,7 +393,8 @@ def test_qc5(self): recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) for i, expected in test_matrix: - recording2 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) + recording2 = np.asarray( + [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) result = hergqc.qc5(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -494,7 +422,8 @@ def test_qc5_1(self): recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) for i, expected in test_matrix: - recording2 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) + recording2 = np.asarray( + [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) result = hergqc.qc5_1(recording1, recording2) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -520,7 +449,8 @@ def test_qc6(self): (100, True), # valc - val = -10.1 ] for i, expected in test_matrix: - recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [0.1 * i] * 500) + recording = np.asarray( + [0, 0.1] * (NOISE_LEN // 2) + [0.1 * i] * 500) result = hergqc.qc6(recording, win=[NOISE_LEN, -1]) self.assertEqual(result[0], expected, f"({i}: {result[1]})") @@ -560,7 +490,7 @@ def test_run_qc(self): after=after, qc_vals_before=qc_vals_before, qc_vals_after=qc_vals_after, - n_sweeps=2, + n_sweeps=self.n_sweeps, ) logging.debug(well, QC.all_passed()) From fc70547e990a763afc150c9e718170734d0a96e2 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Wed, 27 Nov 2024 23:02:24 +0000 Subject: [PATCH 24/33] #21 Run qc3 on all test data --- tests/test_herg_qc.py | 187 +++++++++++++++++++++++++----------------- 1 file changed, 114 insertions(+), 73 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index c9cf864..75e7511 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -56,8 +56,7 @@ def setUp(self): ] plot_dir = "test_output" - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) + os.makedirs(plot_dir, exist_ok=True) self.hergqc = hERGQC( sampling_rate=sampling_rate, @@ -65,17 +64,19 @@ def setUp(self): voltage=self.voltage, ) + def clone_herg_qc(self, plot_dir): + hergqc = copy.deepcopy(self.hergqc) + plot_dir = os.path.join(hergqc.plot_dir, plot_dir) + os.makedirs(plot_dir, exist_ok=True) + hergqc.plot_dir = plot_dir + return hergqc + def test_qc_inputs(self): self.assertTrue(np.all(np.isfinite(self.voltage))) self.assertTrue(np.all(np.isfinite(self.times))) def test_qc1(self): - - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc1") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc1") # qc1 checks that rseal, cm, rseries are within range rseal_lo, rseal_hi = 1e8, 1e12 @@ -183,27 +184,25 @@ def test_qc1(self): ) def test_qc2(self): - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc2") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc2") # qc2 checks that raw and subtracted SNR are above a minimum threshold test_matrix = [ - (10, 8082.1, True), - (1, 74.0, True), - (0.601, 25.1, True), - (0.6, 25.0, False), - (0.5, 16.8, False), - (0.1, 0.5, False), + (10, True, 8082.1), + (1, True, 74.0), + (0.601, True, 25.1), + (0.6, False, 25.0), + (0.599, False, 24.9), + (0.5, False, 16.8), + (0.1, False, 0.5), ] - for i, ex_snr, ex_pass in test_matrix: + for i, ex_pass, ex_snr in test_matrix: recording = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [i] * 500) pass_, snr = hergqc.qc2(recording) self.assertAlmostEqual( - snr, ex_snr, 1, f"QC2: ({i}) {snr} != {ex_snr}") + snr, ex_snr, 1, f"QC2: ({i}) {snr} != {ex_snr}" + ) self.assertEqual(pass_, ex_pass, f"QC2: ({i}) {pass_} != {ex_pass}") # Test on data @@ -237,56 +236,115 @@ def test_qc2(self): ) def test_qc3(self): - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc3") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc3") # qc3 checks that rmsd of two sweeps are similar # Test with same noise, different signal test_matrix = [ - (-10, False), - (-9, False), - (-8, False), # rmsdc - rmsd = -0.6761186497920804 - (-7, True), # rmsdc - rmsd = 0.25355095037585507 - (0, True), # rmsdc - rmsd = 6.761238263598085 - (8, True), # rmsdc - rmsd = 0.6761272774054383 - (9, False), # rmsdc - rmsd = -0.08451158778363332 - (10, False), + (-10, False, -2.5), + (-9, False, -1.6), + (-8, False, -0.7), + (-7, True, 0.3), + (0, True, 6.8), + (8, True, 0.68), + (9, False, -0.08), + (10, False, -0.8), ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) - for i, expected in test_matrix: + for i, ex_pass, ex_rmsd_diff in test_matrix: recording2 = np.asarray( - [0, 0.1] * (NOISE_LEN // 2) + [40 + i] * 500) - result = hergqc.qc3(recording1, recording2) - self.assertEqual(result[0], expected, f"({i}: {result[1]})") + [0, 0.1] * (NOISE_LEN // 2) + [40 + i] * 500 + ) + pass_, rmsd_diff = hergqc.qc3(recording1, recording2) + self.assertAlmostEqual( + rmsd_diff, + ex_rmsd_diff, + 1, + f"QC3: ({i}) {rmsd_diff} != {ex_rmsd_diff}", + ) + self.assertEqual(pass_, ex_pass, f"QC3: ({i}) {pass_} != {ex_pass}") # Test with same signal, different noise - # TODO: Find failing example test_matrix = [ - (10, True), - (100, True), - (1000, True), + (-100, True, 11.3), + (-10, True, 6.3), + (-1, True, 6.7), + (0, True, 6.7), + (1, True, 6.8), + (10, True, 6.4), + (100, True, 11.4), ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) - for i, expected in test_matrix: + for i, ex_pass, ex_rmsd_diff in test_matrix: recording2 = np.asarray( - [0, 0.1 * i] * (NOISE_LEN // 2) + [40] * 500) - result = hergqc.qc3(recording1, recording2) - self.assertEqual(result[0], expected, f"({i}: {result[1]})") + [0, 0.1 * i] * (NOISE_LEN // 2) + [40] * 500 + ) + pass_, rmsd_diff = hergqc.qc3(recording1, recording2) + self.assertAlmostEqual( + rmsd_diff, + ex_rmsd_diff, + 1, + f"QC3: ({i}) {rmsd_diff} != {ex_rmsd_diff}", + ) + self.assertEqual(pass_, ex_pass, f"QC3: ({i}) {pass_} != {ex_pass}") - # TODO: Test on select data + # Test on data + failed_wells_raw = [ + 'A21', 'B05', 'B10', 'C19', 'E09', 'E19', 'F22', 'F23', 'I06', 'K23', + 'L09', 'M05', 'M06', 'M10', 'N12', 'N17', 'O13', 'O15', 'P11' + ] + + failed_wells_E4031 = ['A05', 'A07', 'C19', 'E19', 'J16'] + + failed_wells_subtracted = [ + 'A05', 'A20', 'A21', 'A24', 'B05', 'B07', 'B10', 'B15', 'B21', 'B23', + 'C04', 'C12', 'C14', 'C17', 'C18', 'C19', 'C20', 'D21', 'E04', 'E09', + 'E10', 'E15', 'E16', 'E17', 'E18', 'E19', 'E23', 'F04', 'F06', 'F07', + 'F12', 'F20', 'F21', 'G08', 'G09', 'G10', 'G12', 'G13', 'G16', 'G20', + 'G23', 'G24', 'H03', 'H07', 'H15', 'H19', 'H21', 'H24', 'I04', 'I05', + 'I12', 'I17', 'I21', 'J07', 'J16', 'J19', 'K02', 'K06', 'K22', 'K23', + 'L01', 'L05', 'L08', 'L09', 'L10', 'L11', 'L13', 'L24', 'M01', 'M02', + 'M04', 'M05', 'M06', 'M10', 'M12', 'M21', 'N04', 'N06', 'N08', 'N11', + 'N12', 'N17', 'N20', 'N24', 'O01', 'O03', 'O07', 'O10', 'O13', 'O15', + 'O19', 'O22', 'P01', 'P03', 'P06', 'P08', 'P12', 'P15', 'P17', 'P18' + ] + + for well in self.all_wells: + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) + + pass_raw, rmsd_diff_raw = hergqc.qc3(before[0, :], before[1, :]) + ex_pass_raw = well not in failed_wells_raw + self.assertEqual( + pass_raw, + ex_pass_raw, + f"QC3: {well} (raw) {rmsd_diff_raw}", + ) + + pass_E4031, rmsd_diff_E4031 = hergqc.qc3(after[0, :], after[1, :]) + ex_pass_E4031 = well not in failed_wells_E4031 + self.assertEqual( + pass_E4031, + ex_pass_E4031, + f"QC3: {well} (E4031) {rmsd_diff_E4031}", + ) + + pass_subtracted, rmsd_diff_subtracted = hergqc.qc3( + before[0, :] - after[0, :], + before[1, :] - after[1, :], + ) + ex_pass_subtracted = well not in failed_wells_subtracted + self.assertEqual( + pass_subtracted, + ex_pass_subtracted, + f"QC3: {well} (subtracted) {rmsd_diff_subtracted}", + ) def test_qc4(self): - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc4") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc4") # qc4 checks that rseal, cm, rseries are similar before/after E-4031 change r_lo, r_hi = 1e6, 3e7 @@ -373,11 +431,7 @@ def test_qc4(self): # TODO: Test on select data def test_qc5(self): - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc5") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc5") # qc5 checks that the maximum current during the second half of the # staircase changes by at least 75% of the raw trace after E-4031 addition @@ -401,11 +455,7 @@ def test_qc5(self): # TODO: Test on select data def test_qc5_1(self): - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc5_1") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc5_1") # qc5_1 checks that the RMSD to zero of staircase protocol changes # by at least 50% of the raw trace after E-4031 addition. @@ -430,11 +480,7 @@ def test_qc5_1(self): # TODO: Test on select data def test_qc6(self): - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_qc6") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_qc6") # qc6 checks that the first step up to +40 mV, before the staircase, in # the subtracted trace is bigger than -2 x estimated noise level. @@ -458,12 +504,7 @@ def test_qc6(self): def test_run_qc(self): # Spot check a few wells; could check all, but it's time consuming. - - hergqc = copy.deepcopy(self.hergqc) - plot_dir = os.path.join(hergqc.plot_dir, "test_run_qc") - if not os.path.exists(plot_dir): - os.makedirs(plot_dir) - hergqc.plot_dir = plot_dir + hergqc = self.clone_herg_qc("test_run_qc") test_matrix = [ ("A01", True), From 356469cde15a336231f6dc6dc63afce771ce3caf Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 28 Nov 2024 09:33:33 +0000 Subject: [PATCH 25/33] #21 Run qc4 on all test data --- tests/test_herg_qc.py | 74 ++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 73 insertions(+), 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 75e7511..279bc10 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -428,7 +428,79 @@ def test_qc4(self): f"({i}: {rseals}, {cms}, {rseriess})", ) - # TODO: Test on select data + # Test on data + failed_wells_rseals = [ + 'A04', 'A05', 'A07', 'A16', 'A21', 'A23', 'B02', 'B04', 'B11', 'B16', + 'C10', 'C19', 'C22', 'C23', 'D03', 'D23', 'E01', 'E02', 'E03', 'E07', + 'F23', 'H01', 'H09', 'H17', 'I06', 'I11', 'J11', 'K01', 'K09', 'K12', + 'K14', 'K23', 'M05', 'M10', 'N02', 'N09', 'N17', 'O08', 'O14', 'P16', + 'P24' + ] + + failed_wells_cms = [ + 'A12', 'A13', 'A16', 'A19', 'A20', 'B07', 'B11', 'B13', 'B15', 'B19', + 'B21', 'B23', 'C02', 'C04', 'C07', 'C11', 'C12', 'C14', 'C18', 'D03', + 'D10', 'D14', 'E03', 'E04', 'E07', 'E09', 'E10', 'E15', 'E16', 'E17', + 'E19', 'E22', 'E23', 'F01', 'F03', 'F04', 'F07', 'F12', 'F14', 'F15', + 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', 'G13', 'G16', + 'G20', 'G23', 'G24', 'H01', 'H03', 'H06', 'H07', 'H09', 'H10', 'H15', + 'H19', 'H21', 'H23', 'H24', 'I04', 'I05', 'I07', 'I10', 'I12', 'I16', + 'I17', 'I21', 'J07', 'J16', 'J19', 'J21', 'K02', 'K16', 'K22', 'K23', + 'L01', 'L02', 'L05', 'L08', 'L10', 'L11', 'L13', 'L16', 'L17', 'L18', + 'M01', 'M04', 'M12', 'M15', 'M19', 'M21', 'N06', 'N08', 'N11', 'N14', + 'N18', 'N19', 'N21', 'N24', 'O01', 'O03', 'O07', 'O10', 'O15', 'O17', + 'O19', 'O22', 'O24', 'P01', 'P06', 'P07', 'P08', 'P12', 'P13', 'P14', + 'P15', 'P16', 'P17', 'P18', 'P21', 'P22' + ] + + failed_wells_rseriess = [ + 'A12', 'A13', 'A16', 'A19', 'A20', 'B07', 'B11', 'B13', 'B15', 'B19', + 'B21', 'B23', 'C02', 'C04', 'C07', 'C11', 'C12', 'C14', 'C18', 'D03', + 'D10', 'D14', 'E03', 'E04', 'E07', 'E09', 'E10', 'E15', 'E16', 'E17', + 'E19', 'E22', 'E23', 'F01', 'F03', 'F04', 'F07', 'F12', 'F14', 'F15', + 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', 'G13', 'G16', + 'G20', 'G23', 'G24', 'H01', 'H03', 'H06', 'H07', 'H09', 'H10', 'H15', + 'H19', 'H21', 'H23', 'H24', 'I04', 'I05', 'I07', 'I10', 'I12', 'I16', + 'I17', 'I21', 'J07', 'J16', 'J19', 'J21', 'K02', 'K16', 'K22', 'K23', + 'L01', 'L02', 'L05', 'L08', 'L10', 'L11', 'L13', 'L16', 'L17', 'L18', + 'M01', 'M04', 'M12', 'M15', 'M19', 'M21', 'N06', 'N08', 'N11', 'N14', + 'N18', 'N19', 'N21', 'N24', 'O01', 'O03', 'O07', 'O10', 'O15', 'O17', + 'O19', 'O22', 'O24', 'P01', 'P06', 'P07', 'P08', 'P12', 'P13', 'P14', + 'P15', 'P16', 'P17', 'P18', 'P21', 'P22' + ] + + for well in self.all_wells: + qc_vals_before = np.array(self.qc_vals_before[well])[0, :] + qc_vals_after = np.array(self.qc_vals_after[well])[0, :] + + 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]] + result = hergqc.qc4(rseals, cms, rseriess) + + pass_rseals, d_rseal = result[0] + ex_pass_rseals = well not in failed_wells_rseals + self.assertEqual( + pass_rseals, + ex_pass_rseals, + f"QC4: {well} (rseals) {d_rseal} {rseals}", + ) + + pass_cms, d_cm = result[1] + ex_pass_cms = well not in failed_wells_cms + self.assertEqual( + pass_cms, + ex_pass_cms, + f"QC4: {well} (cms) {d_cm} {cms}", + ) + + pass_rseriess, d_rseries = result[1] + ex_pass_rseriess = well not in failed_wells_rseriess + self.assertEqual( + pass_rseriess, + ex_pass_rseriess, + f"QC4: {well} (rseriess) {d_rseries} {rseriess}", + ) def test_qc5(self): hergqc = self.clone_herg_qc("test_qc5") From 8e1f008f9ee2be0617c74999129b31cfc7172315 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 28 Nov 2024 09:38:57 +0000 Subject: [PATCH 26/33] #21 Fix test_qc4 rseriess --- tests/test_herg_qc.py | 27 ++++++++++++++------------- 1 file changed, 14 insertions(+), 13 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 279bc10..4dc87ce 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -455,18 +455,19 @@ def test_qc4(self): failed_wells_rseriess = [ 'A12', 'A13', 'A16', 'A19', 'A20', 'B07', 'B11', 'B13', 'B15', 'B19', - 'B21', 'B23', 'C02', 'C04', 'C07', 'C11', 'C12', 'C14', 'C18', 'D03', - 'D10', 'D14', 'E03', 'E04', 'E07', 'E09', 'E10', 'E15', 'E16', 'E17', - 'E19', 'E22', 'E23', 'F01', 'F03', 'F04', 'F07', 'F12', 'F14', 'F15', - 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', 'G13', 'G16', - 'G20', 'G23', 'G24', 'H01', 'H03', 'H06', 'H07', 'H09', 'H10', 'H15', - 'H19', 'H21', 'H23', 'H24', 'I04', 'I05', 'I07', 'I10', 'I12', 'I16', - 'I17', 'I21', 'J07', 'J16', 'J19', 'J21', 'K02', 'K16', 'K22', 'K23', - 'L01', 'L02', 'L05', 'L08', 'L10', 'L11', 'L13', 'L16', 'L17', 'L18', - 'M01', 'M04', 'M12', 'M15', 'M19', 'M21', 'N06', 'N08', 'N11', 'N14', - 'N18', 'N19', 'N21', 'N24', 'O01', 'O03', 'O07', 'O10', 'O15', 'O17', - 'O19', 'O22', 'O24', 'P01', 'P06', 'P07', 'P08', 'P12', 'P13', 'P14', - 'P15', 'P16', 'P17', 'P18', 'P21', 'P22' + 'B21', 'B23', 'C02', 'C04', 'C07', 'C11', 'C12', 'C14', 'C18', 'C22', + 'D03', 'D10', 'D14', 'E01', 'E03', 'E04', 'E07', 'E09', 'E10', 'E15', + 'E16', 'E17', 'E22', 'E23', 'F01', 'F03', 'F04', 'F07', 'F12', 'F14', + 'F15', 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', 'G13', + 'G14', 'G16', 'G20', 'G23', 'G24', 'H01', 'H03', 'H04', 'H07', 'H09', + 'H10', 'H13', 'H15', 'H19', 'H21', 'H23', 'H24', 'I04', 'I05', 'I06', + 'I07', 'I10', 'I12', 'I16', 'I17', 'I21', 'J07', 'J14', 'J16', 'J19', + 'J21', 'K02', 'K16', 'K22', 'K23', 'L01', 'L02', 'L05', 'L08', 'L10', + 'L11', 'L12', 'L13', 'L16', 'L17', 'L18', 'M01', 'M02', 'M04', 'M05', + 'M12', 'M13', 'M15', 'M19', 'M20', 'M21', 'N06', 'N08', 'N11', 'N14', + 'N18', 'N19', 'N21', 'N24', 'O01', 'O03', 'O07', 'O10', 'O14', 'O15', + 'O17', 'O19', 'O22', 'O24', 'P01', 'P06', 'P07', 'P08', 'P12', 'P13', + 'P14', 'P15', 'P16', 'P17', 'P18', 'P21', 'P22' ] for well in self.all_wells: @@ -494,7 +495,7 @@ def test_qc4(self): f"QC4: {well} (cms) {d_cm} {cms}", ) - pass_rseriess, d_rseries = result[1] + pass_rseriess, d_rseries = result[2] ex_pass_rseriess = well not in failed_wells_rseriess self.assertEqual( pass_rseriess, From 2b981dc0e6a85bf7412a91ef79700552c4dfb137 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 28 Nov 2024 10:25:40 +0000 Subject: [PATCH 27/33] #21 Expand qc1 data tests --- tests/test_herg_qc.py | 202 +++++++++++++++++++++++++++++++----------- 1 file changed, 151 insertions(+), 51 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 4dc87ce..fc1f81d 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -118,69 +118,169 @@ def test_qc1(self): f"QC1: {rseal}, {cm}, {rseries}", ) - # Test on data - failed_wells_before = [ - 'A10', 'A12', 'A13', 'A16', 'A19', 'A20', 'A24', 'B05', 'B07', 'B11', - 'B12', 'B13', 'B15', 'B19', 'B21', 'B23', 'C02', 'C04', 'C07', 'C09', - 'C11', 'C12', 'C14', 'C18', 'C19', 'C20', 'C23', 'D02', 'D03', 'D05', - 'D09', 'D10', 'D14', 'D15', 'D16', 'D19', 'D21', 'E04', 'E06', 'E07', - 'E10', 'E11', 'E14', 'E15', 'E16', 'E17', 'E19', 'E22', 'E23', 'F01', - 'F03', 'F04', 'F05', 'F07', 'F09', 'F12', 'F13', 'F14', 'F15', 'F18', - 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', 'G13', 'G14', 'G16', - 'G17', 'G20', 'G21', 'G23', 'G24', 'H01', 'H02', 'H03', 'H06', 'H07', - 'H08', 'H10', 'H11', 'H13', 'H14', 'H15', 'H16', 'H19', 'H21', 'H23', - 'H24', 'I01', 'I04', 'I05', 'I06', 'I07', 'I08', 'I10', 'I11', 'I12', - 'I16', 'I17', 'I21', 'J07', 'J10', 'J16', 'J17', 'J19', 'J21', 'J24', - 'K02', 'K03', 'K06', 'K07', 'K11', 'K12', 'K16', 'K17', 'K20', 'K22', - 'K24', 'L01', 'L02', 'L04', 'L05', 'L08', 'L10', 'L11', 'L13', 'L16', - 'L17', 'L18', 'L24', 'M01', 'M04', 'M06', 'M09', 'M12', 'M14', 'M15', - 'M16', 'M19', 'M20', 'M21', 'N03', 'N04', 'N06', 'N08', 'N11', 'N13', - 'N14', 'N16', 'N18', 'N21', 'N24', 'O01', 'O02', 'O03', 'O05', 'O07', - 'O10', 'O11', 'O17', 'O19', 'O22', 'O24', 'P01', 'P03', 'P06', 'P07', - 'P08', 'P09', 'P12', 'P13', 'P14', 'P15', 'P16', 'P17', 'P18', 'P21', - 'P22', 'P24' + # Test on data - values before + failed_wells_rseal_before = [ + 'A10', 'A12', 'A13', 'A16', 'A19', 'A20', 'B05', 'B07', 'B11', 'B12', + 'B13', 'B15', 'B19', 'B21', 'B23', 'C02', 'C04', 'C07', 'C09', 'C11', + 'C12', 'C14', 'C18', 'C19', 'C20', 'D02', 'D03', 'D05', 'D09', 'D10', + 'D14', 'D15', 'D19', 'D21', 'E04', 'E07', 'E10', 'E11', 'E14', 'E15', + 'E16', 'E17', 'E22', 'E23', 'F01', 'F03', 'F04', 'F05', 'F07', 'F09', + 'F12', 'F13', 'F14', 'F15', 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', + 'G09', 'G12', 'G13', 'G14', 'G16', 'G20', 'G23', 'G24', 'H01', 'H02', + 'H03', 'H06', 'H07', 'H08', 'H10', 'H11', 'H13', 'H14', 'H15', 'H16', + 'H19', 'H21', 'H23', 'H24', 'I01', 'I04', 'I05', 'I07', 'I08', 'I10', + 'I11', 'I12', 'I16', 'I17', 'I21', 'J07', 'J10', 'J16', 'J17', 'J19', + 'J21', 'J24', 'K02', 'K03', 'K06', 'K07', 'K11', 'K12', 'K16', 'K17', + 'K20', 'K22', 'K24', 'L01', 'L02', 'L04', 'L05', 'L08', 'L10', 'L11', + 'L13', 'L16', 'L17', 'L18', 'L24', 'M01', 'M04', 'M06', 'M09', 'M12', + 'M16', 'M19', 'M21', 'N03', 'N04', 'N06', 'N08', 'N11', 'N13', 'N14', + 'N16', 'N18', 'N21', 'N24', 'O01', 'O02', 'O03', 'O05', 'O07', 'O10', + 'O11', 'O17', 'O19', 'O22', 'O24', 'P01', 'P03', 'P06', 'P07', 'P08', + 'P09', 'P12', 'P13', 'P14', 'P15', 'P17', 'P18', 'P21', 'P22', 'P24' + ] + + failed_wells_cm_before = [ + 'A12', 'A13', 'A16', 'A19', 'B07', 'B11', 'B13', 'B15', 'B21', 'B23', + 'C02', 'C04', 'C07', 'C11', 'C12', 'C14', 'C18', 'C20', 'D03', 'D10', + 'D14', 'E04', 'E07', 'E10', 'E15', 'E16', 'E17', 'E22', 'E23', 'F01', + 'F03', 'F04', 'F07', 'F12', 'F14', 'F15', 'F18', 'F19', 'F20', 'F21', + 'F24', 'G09', 'G12', 'G13', 'G16', 'G20', 'G23', 'G24', 'H01', 'H03', + 'H06', 'H07', 'H10', 'H15', 'H19', 'H21', 'H23', 'H24', 'I04', 'I05', + 'I07', 'I10', 'I12', 'I16', 'I17', 'I21', 'J07', 'J16', 'J19', 'J21', + 'K02', 'K16', 'K22', 'L01', 'L02', 'L04', 'L05', 'L08', 'L10', 'L11', + 'L13', 'L17', 'L18', 'M01', 'M04', 'M12', 'M19', 'M21', 'N06', 'N08', + 'N11', 'N14', 'N18', 'N21', 'N24', 'O01', 'O03', 'O07', 'O10', 'O17', + 'O19', 'O22', 'O24', 'P01', 'P06', 'P07', 'P08', 'P12', 'P13', 'P14', + 'P15', 'P16', 'P18', 'P21', 'P22' + ] + + failed_wells_rseries_before = [ + 'A12', 'A13', 'A16', 'A19', 'A24', 'B07', 'B11', 'B13', 'B15', 'B21', + 'B23', 'C02', 'C04', 'C07', 'C11', 'C12', 'C14', 'C18', 'C20', 'C23', + 'D03', 'D09', 'D10', 'D14', 'D15', 'D16', 'E04', 'E06', 'E07', 'E10', + 'E15', 'E16', 'E17', 'E19', 'E22', 'E23', 'F01', 'F03', 'F04', 'F07', + 'F12', 'F14', 'F15', 'F18', 'F19', 'F20', 'F21', 'F24', 'G09', 'G12', + 'G13', 'G16', 'G17', 'G20', 'G21', 'G23', 'G24', 'H01', 'H03', 'H07', + 'H10', 'H15', 'H19', 'H21', 'H23', 'H24', 'I04', 'I05', 'I06', 'I07', + 'I10', 'I12', 'I16', 'I17', 'I21', 'J07', 'J16', 'J19', 'J21', 'K02', + 'K16', 'K22', 'L01', 'L02', 'L04', 'L05', 'L08', 'L10', 'L11', 'L13', + 'L17', 'L18', 'M01', 'M04', 'M12', 'M14', 'M15', 'M19', 'M20', 'M21', + 'N06', 'N08', 'N11', 'N14', 'N18', 'N21', 'N24', 'O01', 'O03', 'O07', + 'O10', 'O17', 'O19', 'O22', 'O24', 'P01', 'P06', 'P07', 'P08', 'P12', + 'P13', 'P14', 'P15', 'P16', 'P18', 'P21', 'P22' ] for well in self.all_wells: qc_vals_before = np.array(self.qc_vals_before[well])[0, :] - ex_pass_before = well not in failed_wells_before + result = hergqc.qc1(*qc_vals_before) + + pass_rseal_before, rseal_before = result[0] + ex_pass_rseal_before = well not in failed_wells_rseal_before + self.assertEqual( + pass_rseal_before, + ex_pass_rseal_before, + f"QC1: {well} (rseal before) {rseal_before}", + ) + + pass_cm_before, cm_before = result[1] + ex_pass_cm_before = well not in failed_wells_cm_before + self.assertEqual( + pass_cm_before, + ex_pass_cm_before, + f"QC1: {well} (cm before) {cm_before}", + ) + + pass_rseries_before, rseries_before = result[2] + ex_pass_rseries_before = well not in failed_wells_rseries_before self.assertEqual( - all_passed(hergqc.qc1(*qc_vals_before)), - ex_pass_before, - f"QC1: {well} (before) {qc_vals_before}", + pass_rseries_before, + ex_pass_rseries_before, + f"QC1: {well} (rseries before) {rseries_before}", ) - failed_wells_after = [ - 'A06', 'A08', 'A10', 'A12', 'A13', 'A16', 'A19', 'A20', 'A24', 'B02', - 'B05', 'B07', 'B11', 'B12', 'B13', 'B15', 'B19', 'B21', 'B23', 'C01', - 'C02', 'C04', 'C07', 'C09', 'C11', 'C12', 'C14', 'C18', 'C20', 'C22', - 'C23', 'D03', 'D05', 'D09', 'D10', 'D14', 'D15', 'D19', 'D21', 'E01', - 'E03', 'E04', 'E06', 'E07', 'E09', 'E10', 'E11', 'E14', 'E15', 'E16', - 'E17', 'E19', 'E20', 'E22', 'E23', 'E24', 'F01', 'F03', 'F04', 'F07', - 'F09', 'F12', 'F13', 'F14', 'F15', 'F16', 'F18', 'F19', 'F20', 'F21', - 'F24', 'G06', 'G08', 'G09', 'G12', 'G13', 'G14', 'G15', 'G16', 'G17', - 'G20', 'G21', 'G23', 'G24', 'H01', 'H02', 'H03', 'H04', 'H06', 'H07', - 'H08', 'H09', 'H10', 'H11', 'H13', 'H14', 'H15', 'H16', 'H17', 'H18', - 'H19', 'H20', 'H21', 'H23', 'H24', 'I01', 'I04', 'I05', 'I06', 'I07', - 'I08', 'I10', 'I12', 'I16', 'I17', 'I20', 'I21', 'J03', 'J07', 'J10', - 'J14', 'J16', 'J17', 'J19', 'J21', 'J23', 'J24', 'K02', 'K03', 'K06', - 'K07', 'K14', 'K16', 'K20', 'K22', 'K23', 'K24', 'L01', 'L02', 'L03', - 'L04', 'L05', 'L08', 'L10', 'L11', 'L12', 'L13', 'L16', 'L17', 'L18', - 'L24', 'M01', 'M02', 'M04', 'M05', 'M06', 'M08', 'M09', 'M12', 'M13', - 'M14', 'M15', 'M19', 'M20', 'M21', 'N04', 'N06', 'N08', 'N11', 'N13', - 'N14', 'N16', 'N18', 'N19', 'N21', 'N24', 'O01', 'O02', 'O03', 'O05', - 'O07', 'O08', 'O10', 'O11', 'O14', 'O15', 'O17', 'O19', 'O22', 'O24', - 'P01', 'P03', 'P06', 'P07', 'P08', 'P09', 'P11', 'P12', 'P13', 'P14', - 'P15', 'P16', 'P17', 'P18', 'P19', 'P21', 'P22' + # Test on data - values after + failed_wells_rseal_after = [ + 'A10', 'A12', 'A13', 'A16', 'A19', 'A20', 'A24', 'B02', 'B05', 'B07', + 'B11', 'B12', 'B13', 'B15', 'B21', 'B23', 'C02', 'C04', 'C07', 'C09', + 'C11', 'C12', 'C14', 'C18', 'C20', 'C22', 'D03', 'D05', 'D10', 'D14', + 'D19', 'D21', 'E04', 'E07', 'E10', 'E11', 'E14', 'E15', 'E16', 'E17', + 'E22', 'E23', 'F01', 'F03', 'F04', 'F07', 'F09', 'F12', 'F13', 'F14', + 'F15', 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', 'G08', 'G09', 'G12', + 'G13', 'G15', 'G16', 'G20', 'G23', 'G24', 'H01', 'H03', 'H06', 'H07', + 'H08', 'H09', 'H10', 'H11', 'H14', 'H15', 'H16', 'H17', 'H18', 'H19', + 'H21', 'H23', 'H24', 'I04', 'I05', 'I06', 'I07', 'I08', 'I10', 'I12', + 'I16', 'I17', 'I21', 'J07', 'J10', 'J16', 'J17', 'J19', 'J21', 'J23', + 'J24', 'K02', 'K03', 'K07', 'K14', 'K16', 'K20', 'K22', 'K23', 'K24', + 'L01', 'L02', 'L04', 'L05', 'L08', 'L10', 'L11', 'L13', 'L17', 'L18', + 'L24', 'M01', 'M04', 'M06', 'M09', 'M12', 'M19', 'M21', 'N04', 'N06', + 'N08', 'N11', 'N13', 'N14', 'N16', 'N18', 'N21', 'N24', 'O01', 'O02', + 'O03', 'O07', 'O08', 'O10', 'O11', 'O17', 'O19', 'O22', 'O24', 'P01', + 'P06', 'P07', 'P08', 'P09', 'P12', 'P13', 'P14', 'P15', 'P17', 'P18', + 'P21', 'P22' + ] + + failed_wells_cm_after = [ + 'A12', 'A13', 'A19', 'A20', 'B07', 'B11', 'B13', 'B15', 'B19', 'B21', + 'B23', 'C02', 'C04', 'C11', 'C12', 'C14', 'C18', 'C20', 'D10', 'D14', + 'E03', 'E09', 'E10', 'E15', 'E16', 'E17', 'E19', 'E22', 'E23', 'F01', + 'F03', 'F04', 'F07', 'F12', 'F14', 'F15', 'F18', 'F19', 'F20', 'F21', + 'F24', 'G06', 'G09', 'G12', 'G13', 'G16', 'G20', 'G23', 'G24', 'H01', + 'H03', 'H07', 'H09', 'H10', 'H15', 'H19', 'H21', 'H23', 'H24', 'I04', + 'I05', 'I07', 'I10', 'I12', 'I16', 'I17', 'I21', 'J07', 'J16', 'J19', + 'J21', 'K02', 'K16', 'K22', 'K23', 'L01', 'L02', 'L04', 'L05', 'L08', + 'L10', 'L11', 'L13', 'L16', 'L17', 'L18', 'M01', 'M04', 'M12', 'M15', + 'M19', 'M21', 'N06', 'N08', 'N11', 'N14', 'N18', 'N19', 'N21', 'N24', + 'O01', 'O03', 'O07', 'O10', 'O15', 'O17', 'O19', 'O22', 'O24', 'P01', + 'P06', 'P08', 'P12', 'P13', 'P14', 'P15', 'P16', 'P17', 'P18', 'P21', + 'P22' + ] + + failed_wells_rseries_after = [ + 'A06', 'A08', 'A12', 'A13', 'A19', 'A20', 'A24', 'B07', 'B11', 'B13', + 'B15', 'B19', 'B21', 'B23', 'C01', 'C02', 'C04', 'C09', 'C11', 'C12', + 'C14', 'C18', 'C20', 'C22', 'C23', 'D09', 'D10', 'D14', 'D15', 'D19', + 'E01', 'E03', 'E04', 'E06', 'E09', 'E10', 'E15', 'E16', 'E17', 'E19', + 'E20', 'E22', 'E23', 'E24', 'F01', 'F03', 'F04', 'F07', 'F12', 'F14', + 'F15', 'F16', 'F18', 'F19', 'F20', 'F21', 'F24', 'G06', 'G09', 'G12', + 'G13', 'G14', 'G16', 'G17', 'G20', 'G21', 'G23', 'G24', 'H01', 'H02', + 'H03', 'H04', 'H07', 'H09', 'H10', 'H13', 'H14', 'H15', 'H19', 'H20', + 'H21', 'H23', 'H24', 'I01', 'I04', 'I05', 'I07', 'I10', 'I12', 'I16', + 'I17', 'I20', 'I21', 'J03', 'J07', 'J14', 'J16', 'J19', 'J21', 'K02', + 'K06', 'K16', 'K22', 'K23', 'K24', 'L01', 'L02', 'L03', 'L04', 'L05', + 'L08', 'L10', 'L11', 'L12', 'L13', 'L16', 'L17', 'L18', 'M01', 'M02', + 'M04', 'M05', 'M08', 'M09', 'M12', 'M13', 'M14', 'M15', 'M19', 'M20', + 'M21', 'N06', 'N08', 'N11', 'N14', 'N18', 'N19', 'N21', 'N24', 'O01', + 'O03', 'O05', 'O07', 'O10', 'O11', 'O14', 'O15', 'O17', 'O19', 'O22', + 'O24', 'P01', 'P03', 'P06', 'P08', 'P11', 'P12', 'P13', 'P14', 'P15', + 'P16', 'P17', 'P18', 'P19', 'P21', 'P22' ] for well in self.all_wells: qc_vals_after = np.array(self.qc_vals_after[well])[0, :] - ex_pass_after = well not in failed_wells_after + result = hergqc.qc1(*qc_vals_after) + + pass_rseal_after, rseal_after = result[0] + ex_pass_rseal_after = well not in failed_wells_rseal_after + self.assertEqual( + pass_rseal_after, + ex_pass_rseal_after, + f"QC1: {well} (rseal after) {rseal_after}", + ) + + pass_cm_after, cm_after = result[1] + ex_pass_cm_after = well not in failed_wells_cm_after + self.assertEqual( + pass_cm_after, + ex_pass_cm_after, + f"QC1: {well} (cm after) {cm_after}", + ) + + pass_rseries_after, rseries_after = result[2] + ex_pass_rseries_after = well not in failed_wells_rseries_after self.assertEqual( - all_passed(hergqc.qc1(*qc_vals_after)), - ex_pass_after, - f"QC1: {well} (after) {qc_vals_after}", + pass_rseries_after, + ex_pass_rseries_after, + f"QC1: {well} (rseries after) {rseries_after}", ) def test_qc2(self): From 6bf42debae84fad903157e835315e6453a571589 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 28 Nov 2024 15:19:26 +0000 Subject: [PATCH 28/33] #21 Run qc5 on all test data --- tests/test_herg_qc.py | 138 ++++++++++++++++++++++++++++-------------- 1 file changed, 93 insertions(+), 45 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index fc1f81d..d2c92aa 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -111,10 +111,10 @@ def test_qc1(self): [(0, 0, 0), False], ] - for (rseal, cm, rseries), expected in test_matrix: + for (rseal, cm, rseries), ex_pass in test_matrix: self.assertEqual( all_passed(hergqc.qc1(rseal, cm, rseries)), - expected, + ex_pass, f"QC1: {rseal}, {cm}, {rseries}", ) @@ -353,16 +353,16 @@ def test_qc3(self): ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) - for i, ex_pass, ex_rmsd_diff in test_matrix: + for i, ex_pass, ex_d_rmsd in test_matrix: recording2 = np.asarray( [0, 0.1] * (NOISE_LEN // 2) + [40 + i] * 500 ) - pass_, rmsd_diff = hergqc.qc3(recording1, recording2) + pass_, d_rmsd = hergqc.qc3(recording1, recording2) self.assertAlmostEqual( - rmsd_diff, - ex_rmsd_diff, + d_rmsd, + ex_d_rmsd, 1, - f"QC3: ({i}) {rmsd_diff} != {ex_rmsd_diff}", + f"QC3: ({i}) {d_rmsd} != {ex_d_rmsd}", ) self.assertEqual(pass_, ex_pass, f"QC3: ({i}) {pass_} != {ex_pass}") @@ -378,16 +378,16 @@ def test_qc3(self): ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [40] * 500) - for i, ex_pass, ex_rmsd_diff in test_matrix: + for i, ex_pass, ex_d_rmsd in test_matrix: recording2 = np.asarray( [0, 0.1 * i] * (NOISE_LEN // 2) + [40] * 500 ) - pass_, rmsd_diff = hergqc.qc3(recording1, recording2) + pass_, d_rmsd = hergqc.qc3(recording1, recording2) self.assertAlmostEqual( - rmsd_diff, - ex_rmsd_diff, + d_rmsd, + ex_d_rmsd, 1, - f"QC3: ({i}) {rmsd_diff} != {ex_rmsd_diff}", + f"QC3: ({i}) {d_rmsd} != {ex_d_rmsd}", ) self.assertEqual(pass_, ex_pass, f"QC3: ({i}) {pass_} != {ex_pass}") @@ -416,23 +416,23 @@ def test_qc3(self): before = np.array(self.trace_sweeps_before[well]) after = np.array(self.trace_sweeps_after[well]) - pass_raw, rmsd_diff_raw = hergqc.qc3(before[0, :], before[1, :]) + pass_raw, d_rmsd_raw = hergqc.qc3(before[0, :], before[1, :]) ex_pass_raw = well not in failed_wells_raw self.assertEqual( pass_raw, ex_pass_raw, - f"QC3: {well} (raw) {rmsd_diff_raw}", + f"QC3: {well} (raw) {d_rmsd_raw}", ) - pass_E4031, rmsd_diff_E4031 = hergqc.qc3(after[0, :], after[1, :]) + pass_E4031, d_rmsd_E4031 = hergqc.qc3(after[0, :], after[1, :]) ex_pass_E4031 = well not in failed_wells_E4031 self.assertEqual( pass_E4031, ex_pass_E4031, - f"QC3: {well} (E4031) {rmsd_diff_E4031}", + f"QC3: {well} (E4031) {d_rmsd_E4031}", ) - pass_subtracted, rmsd_diff_subtracted = hergqc.qc3( + pass_subtracted, d_rmsd_subtracted = hergqc.qc3( before[0, :] - after[0, :], before[1, :] - after[1, :], ) @@ -440,7 +440,7 @@ def test_qc3(self): self.assertEqual( pass_subtracted, ex_pass_subtracted, - f"QC3: {well} (subtracted) {rmsd_diff_subtracted}", + f"QC3: {well} (subtracted) {d_rmsd_subtracted}", ) def test_qc4(self): @@ -461,18 +461,18 @@ def test_qc4(self): (5.0, False), ] - for i, expected in test_matrix: + for i, ex_pass in test_matrix: rseals = [r_lo, i * r_lo] self.assertEqual( all_passed(hergqc.qc4(rseals, cms, rseriess)), - expected, + ex_pass, f"({i}: {rseals}, {cms}, {rseriess})", ) rseals = [r_hi, i * r_hi] self.assertEqual( all_passed(hergqc.qc4(rseals, cms, rseriess)), - expected, + ex_pass, f"({i}: {rseals}, {cms}, {rseriess})", ) @@ -487,18 +487,18 @@ def test_qc4(self): (5.0, False), ] - for i, expected in test_matrix: + for i, ex_pass in test_matrix: cms = [c_lo, i * c_lo] self.assertEqual( all_passed(hergqc.qc4(rseals, cms, rseriess)), - expected, + ex_pass, f"({i}: {rseals}, {cms}, {rseriess})", ) cms = [c_hi, i * c_hi] self.assertEqual( all_passed(hergqc.qc4(rseals, cms, rseriess)), - expected, + ex_pass, f"({i}: {rseals}, {cms}, {rseriess})", ) @@ -513,18 +513,18 @@ def test_qc4(self): (5.0, False), ] - for i, expected in test_matrix: + for i, ex_pass in test_matrix: rseriess = [r_lo, i * r_lo] self.assertEqual( all_passed(hergqc.qc4(rseals, cms, rseriess)), - expected, + ex_pass, f"({i}: {rseals}, {cms}, {rseriess})", ) rseriess = [r_hi, i * r_hi] self.assertEqual( all_passed(hergqc.qc4(rseals, cms, rseriess)), - expected, + ex_pass, f"({i}: {rseals}, {cms}, {rseriess})", ) @@ -609,23 +609,71 @@ def test_qc5(self): # qc5 checks that the maximum current during the second half of the # staircase changes by at least 75% of the raw trace after E-4031 addition test_matrix = [ - (-1.0, True), - (0.1, True), - (0.2, True), # max_diffc - max_diff = -0.5 - (0.25, True), # max_diffc - max_diff = 0.0 - (0.3, False), # max_diffc - max_diff = 0.5 - (0.5, False), - (1.0, False), + (-1.0, True, -12.5), + (0.1, True, -1.5), + (0.2, True, -0.5), + (0.24, True, -0.1), + (0.25, True, 0), + (0.26, False, 0.1), + (0.3, False, 0.5), + (0.5, False, 2.5), + (1.0, False, 7.5), ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) - for i, expected in test_matrix: + for i, ex_pass, ex_d_max_diff in test_matrix: recording2 = np.asarray( - [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) - result = hergqc.qc5(recording1, recording2) - self.assertEqual(result[0], expected, f"({i}: {result[1]})") + [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500 + ) + pass_, d_max_diff = hergqc.qc5(recording1, recording2) + self.assertAlmostEqual( + d_max_diff, + ex_d_max_diff, + 1, + f"QC5: ({i}) {d_max_diff} != {ex_d_max_diff}", + ) + self.assertEqual(pass_, ex_pass, f"QC5: ({i}) {pass_} != {ex_pass}") - # TODO: Test on select data + # Test on data + failed_wells = [ + 'A10', 'A12', 'A13', 'A15', 'A19', 'A20', 'A24', 'B05', 'B07', 'B09', + 'B11', 'B12', 'B13', 'B15', 'B18', 'B19', 'B21', 'B23', 'C02', 'C04', + 'C05', 'C07', 'C08', 'C09', 'C11', 'C12', 'C14', 'C17', 'C18', 'C19', + 'C20', 'C21', 'C22', 'C24', 'D02', 'D05', 'D09', 'D10', 'D11', 'D12', + 'D14', 'D17', 'D18', 'D19', 'D21', 'E04', 'E10', 'E11', 'E13', 'E14', + 'E15', 'E16', 'E17', 'E18', 'E19', 'E22', 'E23', 'F01', 'F03', 'F04', + 'F06', 'F07', 'F09', 'F10', 'F12', 'F13', 'F14', 'F15', 'F16', 'F18', + 'F19', 'F20', 'F21', 'F22', 'F24', 'G03', 'G06', 'G08', 'G09', 'G12', + 'G13', 'G14', 'G15', 'G16', 'G18', 'G19', 'G20', 'G23', 'G24', 'H01', + 'H02', 'H03', 'H04', 'H06', 'H07', 'H08', 'H10', 'H11', 'H14', 'H15', + 'H16', 'H18', 'H19', 'H21', 'H23', 'H24', 'I01', 'I04', 'I05', 'I06', + 'I07', 'I08', 'I10', 'I11', 'I12', 'I13', 'I14', 'I16', 'I17', 'I18', + 'I21', 'J07', 'J10', 'J12', 'J15', 'J16', 'J17', 'J18', 'J19', 'J20', + 'J21', 'J23', 'J24', 'K02', 'K03', 'K06', 'K07', 'K11', 'K12', 'K13', + 'K14', 'K16', 'K18', 'K20', 'K22', 'K23', 'K24', 'L01', 'L02', 'L03', + 'L04', 'L05', 'L07', 'L08', 'L10', 'L11', 'L12', 'L13', 'L16', 'L17', + 'L18', 'L24', 'M01', 'M02', 'M03', 'M04', 'M06', 'M08', 'M09', 'M11', + 'M12', 'M14', 'M15', 'M16', 'M17', 'M18', 'M19', 'M21', 'N04', 'N06', + 'N07', 'N08', 'N11', 'N13', 'N14', 'N16', 'N18', 'N20', 'N21', 'N22', + 'N24', 'O01', 'O02', 'O03', 'O07', 'O08', 'O10', 'O12', 'O15', 'O16', + 'O17', 'O18', 'O19', 'O20', 'O21', 'O22', 'O24', 'P01', 'P03', 'P05', + 'P06', 'P07', 'P08', 'P09', 'P10', 'P12', 'P13', 'P14', 'P15', 'P17', + 'P18', 'P20', 'P21', 'P22', 'P24' + ] + + for well in self.all_wells: + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) + + pass_, d_max_diff = hergqc.qc5( + before[0, :], after[0, :], hergqc.qc5_win + ) + ex_pass = well not in failed_wells + self.assertEqual( + pass_, + ex_pass, + f"QC5: {well} {d_max_diff}", + ) def test_qc5_1(self): hergqc = self.clone_herg_qc("test_qc5_1") @@ -644,11 +692,11 @@ def test_qc5_1(self): ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) - for i, expected in test_matrix: + for i, ex_pass in test_matrix: recording2 = np.asarray( [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) result = hergqc.qc5_1(recording1, recording2) - self.assertEqual(result[0], expected, f"({i}: {result[1]})") + self.assertEqual(result[0], ex_pass, f"({i}: {result[1]})") # TODO: Test on select data @@ -667,11 +715,11 @@ def test_qc6(self): (10, True), # valc - val = -1.1 (100, True), # valc - val = -10.1 ] - for i, expected in test_matrix: + for i, ex_pass in test_matrix: recording = np.asarray( [0, 0.1] * (NOISE_LEN // 2) + [0.1 * i] * 500) result = hergqc.qc6(recording, win=[NOISE_LEN, -1]) - self.assertEqual(result[0], expected, f"({i}: {result[1]})") + self.assertEqual(result[0], ex_pass, f"({i}: {result[1]})") # TODO: Test on select data @@ -688,7 +736,7 @@ def test_run_qc(self): ("D01", False), ] - for well, expected in test_matrix: + for well, ex_pass in test_matrix: with self.subTest(well): # Take values from the first sweep only before = np.array(self.trace_sweeps_before[well]) @@ -714,4 +762,4 @@ def test_run_qc(self): if not QC.qc_passed(label): trace += f"{well} {label}: {result}\n" - self.assertEqual(QC.all_passed(), expected, trace) + self.assertEqual(QC.all_passed(), ex_pass, trace) From 467fae21e4a83122ef75d2fb6385eee3286b637d Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Thu, 28 Nov 2024 16:58:33 +0000 Subject: [PATCH 29/33] #21 Run qc5_1 on all test data --- tests/test_herg_qc.py | 77 +++++++++++++++++++++++++++++++++++-------- 1 file changed, 64 insertions(+), 13 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index d2c92aa..345a8b3 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -681,24 +681,75 @@ def test_qc5_1(self): # qc5_1 checks that the RMSD to zero of staircase protocol changes # by at least 50% of the raw trace after E-4031 addition. test_matrix = [ - (-1.0, False), # rmsd0_diffc - rmsd0_diff = 4.2 - (-0.5, False), # rmsd0_diffc - rmsd0_diff = 0.0001 - (-0.4, True), # rmsd0_diffc - rmsd0_diff = -0.8 - (-0.1, True), # rmsd0_diffc - rmsd0_diff = -3.4 - (0.1, True), # rmsd0_diffc - rmsd0_diff = -3.4 - (0.4, True), # rmsd0_diffc - rmsd0_diff = -0.8 - (0.5, False), # rmsd0_diffc - rmsd0_diff = 0.0001 - (1.0, False), # rmsd0_diffc - rmsd0_diff = 4.2 + (-1.0, False, 4.22), + (-0.5, False, 0), + (-0.412, True, -0.74), + (-0.4, True, -0.84), + (-0.1, True, -3.38), + (0.1, True, -3.38), + (0.4, True, -0.84), + (0.412, True, -0.74), + (0.5, False, 0), + (1.0, False, 4.22), ] recording1 = np.asarray([0, 0.1] * (NOISE_LEN // 2) + [10] * 500) - for i, ex_pass in test_matrix: + for i, ex_pass, ex_d_max_diff in test_matrix: recording2 = np.asarray( - [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500) - result = hergqc.qc5_1(recording1, recording2) - self.assertEqual(result[0], ex_pass, f"({i}: {result[1]})") + [0, 0.1] * (NOISE_LEN // 2) + [10 * i] * 500 + ) + pass_, d_max_diff = hergqc.qc5_1(recording1, recording2) + self.assertAlmostEqual( + d_max_diff, + ex_d_max_diff, + 2, + f"QC5_1: ({i}) {d_max_diff} != {ex_d_max_diff}", + ) + self.assertEqual( + pass_, ex_pass, f"QC5_1: ({i}) {pass_} != {ex_pass}" + ) - # TODO: Test on select data + # Test on data + failed_wells = [ + 'A05', 'A10', 'A12', 'A13', 'A15', 'A19', 'A20', 'A24', 'B02', 'B05', + 'B07', 'B09', 'B10', 'B12', 'B13', 'B14', 'B15', 'B18', 'B19', 'B21', + 'B23', 'C02', 'C04', 'C05', 'C07', 'C08', 'C09', 'C11', 'C12', 'C14', + 'C17', 'C18', 'C20', 'C21', 'C22', 'C24', 'D02', 'D04', 'D05', 'D09', + 'D10', 'D11', 'D12', 'D13', 'D16', 'D17', 'D18', 'D19', 'D21', 'E04', + 'E09', 'E10', 'E11', 'E13', 'E14', 'E15', 'E16', 'E17', 'E18', 'E19', + 'E21', 'E22', 'E23', 'F01', 'F02', 'F03', 'F04', 'F05', 'F06', 'F07', + 'F09', 'F10', 'F11', 'F12', 'F13', 'F14', 'F15', 'F16', 'F18', 'F19', + 'F20', 'F21', 'F22', 'F24', 'G03', 'G06', 'G08', 'G09', 'G10', 'G12', + 'G13', 'G14', 'G15', 'G16', 'G18', 'G19', 'G20', 'G23', 'G24', 'H01', + 'H02', 'H03', 'H04', 'H06', 'H07', 'H08', 'H09', 'H10', 'H11', 'H13', + 'H14', 'H15', 'H16', 'H17', 'H18', 'H19', 'H21', 'H23', 'H24', 'I01', + 'I03', 'I04', 'I05', 'I06', 'I07', 'I08', 'I10', 'I12', 'I13', 'I14', + 'I15', 'I16', 'I17', 'I18', 'I19', 'I20', 'I21', 'J06', 'J07', 'J08', + 'J10', 'J12', 'J15', 'J16', 'J17', 'J18', 'J19', 'J20', 'J21', 'J23', + 'J24', 'K01', 'K02', 'K03', 'K05', 'K06', 'K07', 'K10', 'K11', 'K13', + 'K14', 'K16', 'K18', 'K20', 'K22', 'K23', 'K24', 'L01', 'L02', 'L03', + 'L04', 'L05', 'L06', 'L07', 'L08', 'L10', 'L11', 'L12', 'L13', 'L16', + 'L17', 'L18', 'L20', 'L24', 'M01', 'M02', 'M03', 'M04', 'M06', 'M07', + 'M08', 'M09', 'M11', 'M12', 'M14', 'M16', 'M17', 'M18', 'M19', 'M21', + 'N01', 'N03', 'N04', 'N06', 'N07', 'N08', 'N11', 'N13', 'N14', 'N16', + 'N18', 'N20', 'N21', 'N22', 'N23', 'N24', 'O01', 'O02', 'O03', 'O04', + 'O05', 'O06', 'O07', 'O08', 'O10', 'O11', 'O12', 'O15', 'O16', 'O17', + 'O18', 'O19', 'O20', 'O21', 'O22', 'O24', 'P01', 'P03', 'P05', 'P06', + 'P07', 'P08', 'P09', 'P10', 'P12', 'P13', 'P14', 'P15', 'P16', 'P17', + 'P18', 'P20', 'P21', 'P22' + ] + + for well in self.all_wells: + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) + + pass_, d_rmsd = hergqc.qc5_1(before[0, :], after[0, :], label='1') + ex_pass = well not in failed_wells + self.assertEqual( + pass_, + ex_pass, + f"QC5_1: {well} {d_rmsd}", + ) def test_qc6(self): hergqc = self.clone_herg_qc("test_qc6") From f2ac76005d1ed1a9f1a9cf8dcd25ed57c9c033b5 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Fri, 29 Nov 2024 09:57:56 +0000 Subject: [PATCH 30/33] #21 Run qc6 on all test data --- tests/test_herg_qc.py | 230 +++++++++++++++++++++++++++++++++--------- 1 file changed, 183 insertions(+), 47 deletions(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index 345a8b3..d8f6ae8 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -63,6 +63,7 @@ def setUp(self): plot_dir=plot_dir, voltage=self.voltage, ) + self.hergqc._debug = False # Turn this on to output plots def clone_herg_qc(self, plot_dir): hergqc = copy.deepcopy(self.hergqc) @@ -72,8 +73,19 @@ def clone_herg_qc(self, plot_dir): return hergqc def test_qc_inputs(self): - self.assertTrue(np.all(np.isfinite(self.voltage))) - self.assertTrue(np.all(np.isfinite(self.times))) + times = self.times + voltage = self.hergqc.voltage + qc6_win = self.hergqc.qc6_win + qc6_1_win = self.hergqc.qc6_1_win + qc6_2_win = self.hergqc.qc6_2_win + + self.assertTrue(np.all(np.isfinite(times))) + self.assertTrue(np.all(np.isfinite(voltage))) + + # Ensure thats the windows are correct by checking the voltage trace + assert np.all(np.abs(voltage[qc6_win[0]: qc6_win[1]] - 40.0)) < 1e-8 + assert np.all(np.abs(voltage[qc6_1_win[0]: qc6_1_win[1]] - 40.0)) < 1e-8 + assert np.all(np.abs(voltage[qc6_2_win[0]: qc6_2_win[1]] - 40.0)) < 1e-8 def test_qc1(self): hergqc = self.clone_herg_qc("test_qc1") @@ -757,60 +769,184 @@ def test_qc6(self): # qc6 checks that the first step up to +40 mV, before the staircase, in # the subtracted trace is bigger than -2 x estimated noise level. test_matrix = [ - (-100, False), # valc - val = 9.9 - (-10, False), # valc - val = 0.9 - (-2, False), # valc - val = 0.1 - (-1, True), # valc - val = 0 - (1, True), # valc - val = -0.2 - (2, True), # valc - val = -0.3 - (10, True), # valc - val = -1.1 - (100, True), # valc - val = -10.1 + (-100, False, 9.9), + (-10, False, 0.9), + (-1, True, 0), + (0, True, -0.1), + (1, True, -0.2), + (10, True, -1.1), + (100, True, -10.1), ] - for i, ex_pass in test_matrix: + + for i, ex_pass, ex_d_val in test_matrix: recording = np.asarray( - [0, 0.1] * (NOISE_LEN // 2) + [0.1 * i] * 500) - result = hergqc.qc6(recording, win=[NOISE_LEN, -1]) - self.assertEqual(result[0], ex_pass, f"({i}: {result[1]})") + [0, 0.1] * (NOISE_LEN // 2) + [0.1 * i] * 500 + ) + pass_, d_val = hergqc.qc6(recording, win=[NOISE_LEN, -1]) + self.assertAlmostEqual( + d_val, + ex_d_val, + 1, + f"QC6: ({i}) {d_val} != {ex_d_val}", + ) + self.assertEqual(pass_, ex_pass, f"QC6: ({i}) {pass_} != {ex_pass}") + + # Test on data + failed_wells_0 = [ + 'A05', 'A07', 'A11', 'A12', 'A13', 'A20', 'A22', 'A24', 'B02', 'B05', + 'B07', 'B09', 'B10', 'B15', 'B23', 'C02', 'C04', 'C14', 'C18', 'C22', + 'D01', 'D12', 'E04', 'E09', 'E10', 'E15', 'E16', 'E17', 'E19', 'E21', + 'F02', 'F03', 'F07', 'F10', 'F18', 'F20', 'F21', 'F22', 'G06', 'G08', + 'G10', 'G15', 'G16', 'G20', 'G23', 'H01', 'H09', 'H11', 'H15', 'H17', + 'H19', 'H21', 'H23', 'H24', 'I04', 'I06', 'I10', 'I12', 'I15', 'I21', + 'J07', 'J09', 'J16', 'J21', 'K01', 'K02', 'K03', 'K09', 'K13', 'K14', + 'K18', 'K22', 'K23', 'L01', 'L02', 'L05', 'L07', 'L08', 'L10', 'L11', + 'L13', 'L17', 'L18', 'L21', 'L23', 'L24', 'M01', 'M04', 'M06', 'M07', + 'M11', 'M12', 'M17', 'N06', 'N11', 'N14', 'N20', 'N24', 'O01', 'O07', + 'O08', 'O16', 'O17', 'O19', 'O22', 'O24', 'P01', 'P06', 'P08', 'P12', + 'P14', 'P15', 'P17', 'P18', 'P21' + ] + + failed_wells_1 = [ + 'A05', 'A07', 'A11', 'A12', 'A13', 'A20', 'A22', 'A24', 'B02', 'B05', + 'B07', 'B10', 'B15', 'C02', 'C04', 'C14', 'C18', 'C22', 'D01', 'D12', + 'E04', 'E09', 'E10', 'E14', 'E15', 'E16', 'E19', 'E21', 'F02', 'F03', + 'F07', 'F10', 'F12', 'F18', 'F20', 'F21', 'F22', 'G06', 'G08', 'G09', + 'G10', 'G15', 'G16', 'G20', 'G23', 'H01', 'H09', 'H11', 'H15', 'H17', + 'H19', 'H21', 'H23', 'H24', 'I04', 'I06', 'I10', 'I12', 'I15', 'I21', + 'J07', 'J16', 'J21', 'K01', 'K02', 'K09', 'K13', 'K14', 'K23', 'L01', + 'L02', 'L05', 'L07', 'L08', 'L10', 'L11', 'L13', 'L17', 'L18', 'L21', + 'L24', 'M01', 'M04', 'M06', 'M07', 'M11', 'M12', 'M17', 'N06', 'N11', + 'N14', 'N20', 'N24', 'O01', 'O06', 'O07', 'O08', 'O17', 'O19', 'O22', + 'O24', 'P01', 'P06', 'P08', 'P12', 'P14', 'P15', 'P17', 'P18', 'P21' + ] + + failed_wells_2 = [ + 'A05', 'A07', 'A10', 'A11', 'A12', 'A13', 'A20', 'A22', 'A24', 'B02', + 'B05', 'B07', 'B10', 'B15', 'C02', 'C04', 'C14', 'C22', 'D01', 'D12', + 'E04', 'E09', 'E10', 'E14', 'E15', 'E16', 'E21', 'F02', 'F03', 'F07', + 'F10', 'F12', 'F18', 'F20', 'F21', 'F22', 'G06', 'G08', 'G09', 'G10', + 'G15', 'G16', 'G20', 'G23', 'H01', 'H07', 'H09', 'H11', 'H15', 'H17', + 'H19', 'H21', 'H23', 'H24', 'I06', 'I10', 'I12', 'I15', 'I21', 'J07', + 'J16', 'J21', 'K01', 'K02', 'K09', 'K13', 'K14', 'K23', 'L01', 'L02', + 'L08', 'L10', 'L11', 'L13', 'L17', 'L18', 'L21', 'L24', 'M01', 'M04', + 'M06', 'M07', 'M11', 'M12', 'M17', 'N06', 'N11', 'N14', 'N24', 'O01', + 'O06', 'O07', 'O08', 'O17', 'O19', 'O24', 'P01', 'P06', 'P08', 'P12', + 'P14', 'P15', 'P17', 'P18', 'P21' + ] + + for well in self.all_wells: + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) + + subtracted_0 = [] + subtracted_1 = [] + subtracted_2 = [] + + for i in range(before.shape[0]): + subtracted_0.append( + hergqc.qc6( + (before[i, :] - after[i, :]), hergqc.qc6_win, label="0" + ) + ) + + subtracted_1.append( + hergqc.qc6( + (before[i, :] - after[i, :]), + hergqc.qc6_1_win, + label="1", + ) + ) + + subtracted_2.append( + hergqc.qc6( + (before[i, :] - after[i, :]), + hergqc.qc6_2_win, + label="2", + ) + ) + + ex_pass_0 = well not in failed_wells_0 + self.assertEqual( + all_passed(subtracted_0), + ex_pass_0, + f"QC6 (0): {well} {subtracted_0}", + ) - # TODO: Test on select data + ex_pass_1 = well not in failed_wells_1 + self.assertEqual( + all_passed(subtracted_1), + ex_pass_1, + f"QC6 (1): {well} {subtracted_1}", + ) + + ex_pass_2 = well not in failed_wells_2 + self.assertEqual( + all_passed(subtracted_2), + ex_pass_2, + f"QC6 (2): {well} {subtracted_2}", + ) def test_run_qc(self): - # Spot check a few wells; could check all, but it's time consuming. + # Test all wells hergqc = self.clone_herg_qc("test_run_qc") - test_matrix = [ - ("A01", True), - ("A02", True), - ("A03", True), - ("A04", False), - ("A05", False), - ("D01", False), + failed_wells = [ + 'A04', 'A05', 'A06', 'A07', 'A08', 'A10', 'A11', 'A12', 'A13', 'A15', + 'A16', 'A19', 'A20', 'A21', 'A22', 'A23', 'A24', 'B02', 'B04', 'B05', + 'B07', 'B09', 'B10', 'B11', 'B12', 'B13', 'B14', 'B15', 'B16', 'B18', + 'B19', 'B21', 'B23', 'C01', 'C02', 'C04', 'C05', 'C07', 'C08', 'C09', + 'C10', 'C11', 'C12', 'C14', 'C17', 'C18', 'C19', 'C20', 'C21', 'C22', + 'C23', 'C24', 'D01', 'D02', 'D03', 'D04', 'D05', 'D09', 'D10', 'D11', + 'D12', 'D13', 'D14', 'D15', 'D16', 'D17', 'D18', 'D19', 'D21', 'D23', + 'E01', 'E02', 'E03', 'E04', 'E06', 'E07', 'E09', 'E10', 'E11', 'E13', + 'E14', 'E15', 'E16', 'E17', 'E18', 'E19', 'E20', 'E21', 'E22', 'E23', + 'E24', 'F01', 'F02', 'F03', 'F04', 'F05', 'F06', 'F07', 'F09', 'F10', + 'F11', 'F12', 'F13', 'F14', 'F15', 'F16', 'F18', 'F19', 'F20', 'F21', + 'F22', 'F23', 'F24', 'G03', 'G06', 'G08', 'G09', 'G10', 'G12', 'G13', + 'G14', 'G15', 'G16', 'G17', 'G18', 'G19', 'G20', 'G21', 'G23', 'G24', + 'H01', 'H02', 'H03', 'H04', 'H06', 'H07', 'H08', 'H09', 'H10', 'H11', + 'H13', 'H14', 'H15', 'H16', 'H17', 'H18', 'H19', 'H20', 'H21', 'H23', + 'H24', 'I01', 'I03', 'I04', 'I05', 'I06', 'I07', 'I08', 'I10', 'I11', + 'I12', 'I13', 'I14', 'I15', 'I16', 'I17', 'I18', 'I19', 'I20', 'I21', + 'J03', 'J06', 'J07', 'J08', 'J09', 'J10', 'J11', 'J12', 'J14', 'J15', + 'J16', 'J17', 'J18', 'J19', 'J20', 'J21', 'J23', 'J24', 'K01', 'K02', + 'K03', 'K05', 'K06', 'K07', 'K09', 'K10', 'K11', 'K12', 'K13', 'K14', + 'K16', 'K17', 'K18', 'K20', 'K22', 'K23', 'K24', 'L01', 'L02', 'L03', + 'L04', 'L05', 'L06', 'L07', 'L08', 'L10', 'L11', 'L12', 'L13', 'L16', + 'L17', 'L18', 'L20', 'L21', 'L23', 'L24', 'M01', 'M02', 'M03', 'M04', + 'M05', 'M06', 'M07', 'M08', 'M09', 'M10', 'M11', 'M12', 'M13', 'M14', + 'M15', 'M16', 'M17', 'M18', 'M19', 'M20', 'M21', 'N01', 'N02', 'N03', + 'N04', 'N06', 'N07', 'N08', 'N09', 'N11', 'N13', 'N14', 'N16', 'N17', + 'N18', 'N19', 'N20', 'N21', 'N22', 'N23', 'N24', 'O01', 'O02', 'O03', + 'O04', 'O05', 'O06', 'O07', 'O08', 'O10', 'O11', 'O12', 'O13', 'O14', + 'O15', 'O16', 'O17', 'O18', 'O19', 'O20', 'O21', 'O22', 'O24', 'P01', + 'P03', 'P05', 'P06', 'P07', 'P08', 'P09', 'P10', 'P11', 'P12', 'P13', + 'P14', 'P15', 'P16', 'P17', 'P18', 'P19', 'P20', 'P21', 'P22', 'P24' ] - for well, ex_pass in test_matrix: - with self.subTest(well): - # Take values from the first sweep only - before = np.array(self.trace_sweeps_before[well]) - after = np.array(self.trace_sweeps_after[well]) - - qc_vals_before = np.array(self.qc_vals_before[well])[0, :] - qc_vals_after = np.array(self.qc_vals_after[well])[0, :] - - QC = hergqc.run_qc( - voltage_steps=self.voltage_steps, - times=self.times, - before=before, - after=after, - qc_vals_before=qc_vals_before, - qc_vals_after=qc_vals_after, - n_sweeps=self.n_sweeps, - ) + for well in self.all_wells: + before = np.array(self.trace_sweeps_before[well]) + after = np.array(self.trace_sweeps_after[well]) - logging.debug(well, QC.all_passed()) + # Take values from the first sweep only + qc_vals_before = np.array(self.qc_vals_before[well])[0, :] + qc_vals_after = np.array(self.qc_vals_after[well])[0, :] - trace = "" - for label, result in QC.items(): - if not QC.qc_passed(label): - trace += f"{well} {label}: {result}\n" + QC = hergqc.run_qc( + voltage_steps=self.voltage_steps, + times=self.times, + before=before, + after=after, + qc_vals_before=qc_vals_before, + qc_vals_after=qc_vals_after, + n_sweeps=self.n_sweeps, + ) - self.assertEqual(QC.all_passed(), ex_pass, trace) + trace = "" + for label, result in QC.items(): + if not QC.qc_passed(label): + trace += f"{well} {label}: {result}\n" + + ex_pass = well not in failed_wells + self.assertEqual(QC.all_passed(), ex_pass, trace) From 884346bf7fde403ef34e895603ce76a11a3a8849 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Fri, 29 Nov 2024 10:16:02 +0000 Subject: [PATCH 31/33] #21 Add docstring for QCDict --- pcpostprocess/hergQC.py | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index f07387d..6b576bd 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -10,6 +10,15 @@ class QCDict: + """ + Stores the results from QC checks. + + Each entry is a label -> [(bool, value),...] mapping. + The bool in each tuple indicates whether the QC passed, + and the value is the result that was checked (e.g. the SNR value). + The list can contain multiple tuples if the QC checks multiple values, as + QC1 does, or if the checks are run multiple times e.g. once per sweep. + """ labels = [ "qc1.rseal", @@ -31,7 +40,9 @@ class QCDict: ] def __init__(self): - self._dict = OrderedDict([(label, [(False, None)]) for label in QCDict.labels]) + self._dict = OrderedDict( + [(label, [(False, None)]) for label in QCDict.labels] + ) def __str__(self): return self._dict.__str__() From f40f8f33ae63bd9790fde210d4f7d050a17354b7 Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Fri, 29 Nov 2024 10:18:08 +0000 Subject: [PATCH 32/33] #21 Linter fixes --- tests/test_herg_qc.py | 1 - 1 file changed, 1 deletion(-) diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index d8f6ae8..a831945 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -1,5 +1,4 @@ import copy -import logging import os import string import unittest From 0bbfd445089ae467d61c3c75d30af71801ef349f Mon Sep 17 00:00:00 2001 From: Kwabena N Amponsah Date: Fri, 29 Nov 2024 17:30:05 +0000 Subject: [PATCH 33/33] #21 Implement QCDict with UserDict --- pcpostprocess/hergQC.py | 42 +++++++++++++---------------------------- tests/test_herg_qc.py | 33 ++++++++++++++++---------------- 2 files changed, 29 insertions(+), 46 deletions(-) diff --git a/pcpostprocess/hergQC.py b/pcpostprocess/hergQC.py index 6b576bd..e248a3c 100644 --- a/pcpostprocess/hergQC.py +++ b/pcpostprocess/hergQC.py @@ -1,6 +1,6 @@ import logging import os -from collections import OrderedDict +from collections import OrderedDict, UserDict import matplotlib.pyplot as plt import numpy as np @@ -9,7 +9,7 @@ NOISE_LEN = 200 -class QCDict: +class QCDict(UserDict): """ Stores the results from QC checks. @@ -40,36 +40,16 @@ class QCDict: ] def __init__(self): - self._dict = OrderedDict( - [(label, [(False, None)]) for label in QCDict.labels] - ) - - def __str__(self): - return self._dict.__str__() - - def __repr__(self): - return self._dict.__repr__() - - def __getitem__(self, key): - return self._dict.__getitem__(key) + super().__init__([(label, [(False, None)]) for label in QCDict.labels]) def __setitem__(self, key, value): if key not in QCDict.labels: raise KeyError(f"Invalid QC key: {key}") - self._dict.__setitem__(key, value) - - def keys(self): - return self._dict.keys() - - def items(self): - return self._dict.items() - - def values(self): - return self._dict.values() + super().__setitem__(key, value) def qc_passed(self, label): """Return whether a single QC passed.""" - return all([x for x, _ in self._dict[label]]) + return all([x for x, _ in self[label]]) def passed_list(self): """Return a list of booleans indicating whether each QC passed.""" @@ -250,9 +230,12 @@ def run_qc(self, voltage_steps, times, QC['qc6.1.subtracted'] = [] QC['qc6.2.subtracted'] = [] for i in range(before.shape[0]): - 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") + 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) @@ -463,7 +446,8 @@ 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}") + self.logger.debug( + f"rmsd0_diff: {rmsd0_diff}, rmsd0c: {rmsd0_diffc}") result = False else: result = True diff --git a/tests/test_herg_qc.py b/tests/test_herg_qc.py index a831945..9a06017 100644 --- a/tests/test_herg_qc.py +++ b/tests/test_herg_qc.py @@ -64,7 +64,7 @@ def setUp(self): ) self.hergqc._debug = False # Turn this on to output plots - def clone_herg_qc(self, plot_dir): + def clone_hergqc(self, plot_dir): hergqc = copy.deepcopy(self.hergqc) plot_dir = os.path.join(hergqc.plot_dir, plot_dir) os.makedirs(plot_dir, exist_ok=True) @@ -81,15 +81,15 @@ def test_qc_inputs(self): self.assertTrue(np.all(np.isfinite(times))) self.assertTrue(np.all(np.isfinite(voltage))) - # Ensure thats the windows are correct by checking the voltage trace + # Ensures that the windows are correct by checking the voltage trace assert np.all(np.abs(voltage[qc6_win[0]: qc6_win[1]] - 40.0)) < 1e-8 assert np.all(np.abs(voltage[qc6_1_win[0]: qc6_1_win[1]] - 40.0)) < 1e-8 assert np.all(np.abs(voltage[qc6_2_win[0]: qc6_2_win[1]] - 40.0)) < 1e-8 def test_qc1(self): - hergqc = self.clone_herg_qc("test_qc1") - # qc1 checks that rseal, cm, rseries are within range + hergqc = self.clone_hergqc("test_qc1") + rseal_lo, rseal_hi = 1e8, 1e12 rseal_mid = (rseal_lo + rseal_hi) / 2 rseal_tol = 0.1 @@ -295,9 +295,9 @@ def test_qc1(self): ) def test_qc2(self): - hergqc = self.clone_herg_qc("test_qc2") - # qc2 checks that raw and subtracted SNR are above a minimum threshold + hergqc = self.clone_hergqc("test_qc2") + test_matrix = [ (10, True, 8082.1), (1, True, 74.0), @@ -347,9 +347,8 @@ def test_qc2(self): ) def test_qc3(self): - hergqc = self.clone_herg_qc("test_qc3") - # qc3 checks that rmsd of two sweeps are similar + hergqc = self.clone_hergqc("test_qc3") # Test with same noise, different signal test_matrix = [ @@ -455,9 +454,9 @@ def test_qc3(self): ) def test_qc4(self): - hergqc = self.clone_herg_qc("test_qc4") - # qc4 checks that rseal, cm, rseries are similar before/after E-4031 change + hergqc = self.clone_hergqc("test_qc4") + r_lo, r_hi = 1e6, 3e7 c_lo, c_hi = 1e-12, 1e-10 @@ -615,10 +614,10 @@ def test_qc4(self): ) def test_qc5(self): - hergqc = self.clone_herg_qc("test_qc5") - # qc5 checks that the maximum current during the second half of the # staircase changes by at least 75% of the raw trace after E-4031 addition + hergqc = self.clone_hergqc("test_qc5") + test_matrix = [ (-1.0, True, -12.5), (0.1, True, -1.5), @@ -687,10 +686,10 @@ def test_qc5(self): ) def test_qc5_1(self): - hergqc = self.clone_herg_qc("test_qc5_1") - # qc5_1 checks that the RMSD to zero of staircase protocol changes # by at least 50% of the raw trace after E-4031 addition. + hergqc = self.clone_hergqc("test_qc5_1") + test_matrix = [ (-1.0, False, 4.22), (-0.5, False, 0), @@ -763,10 +762,10 @@ def test_qc5_1(self): ) def test_qc6(self): - hergqc = self.clone_herg_qc("test_qc6") - # qc6 checks that the first step up to +40 mV, before the staircase, in # the subtracted trace is bigger than -2 x estimated noise level. + hergqc = self.clone_hergqc("test_qc6") + test_matrix = [ (-100, False, 9.9), (-10, False, 0.9), @@ -888,7 +887,7 @@ def test_qc6(self): def test_run_qc(self): # Test all wells - hergqc = self.clone_herg_qc("test_run_qc") + hergqc = self.clone_hergqc("test_run_qc") failed_wells = [ 'A04', 'A05', 'A06', 'A07', 'A08', 'A10', 'A11', 'A12', 'A13', 'A15',