Coverage for tasks/cardinal_expdetthreshold.py : 34%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1#!/usr/bin/env python
3"""
4camcops_server/tasks/cardinal_expdetthreshold.py
6===============================================================================
8 Copyright (C) 2012-2020 Rudolf Cardinal (rudolf@pobox.com).
10 This file is part of CamCOPS.
12 CamCOPS is free software: you can redistribute it and/or modify
13 it under the terms of the GNU General Public License as published by
14 the Free Software Foundation, either version 3 of the License, or
15 (at your option) any later version.
17 CamCOPS is distributed in the hope that it will be useful,
18 but WITHOUT ANY WARRANTY; without even the implied warranty of
19 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 GNU General Public License for more details.
22 You should have received a copy of the GNU General Public License
23 along with CamCOPS. If not, see <https://www.gnu.org/licenses/>.
25===============================================================================
27"""
29import math
30import logging
31from typing import List, Optional, Type
33from cardinal_pythonlib.maths_numpy import inv_logistic, logistic
34import cardinal_pythonlib.rnc_web as ws
35import numpy as np
36from sqlalchemy.sql.schema import Column
37from sqlalchemy.sql.sqltypes import Float, Integer, Text, UnicodeText
39from camcops_server.cc_modules.cc_constants import (
40 CssClass,
41 MatplotlibConstants,
42 PlotDefaults,
43)
44from camcops_server.cc_modules.cc_db import (
45 ancillary_relationship,
46 GenericTabletRecordMixin,
47 TaskDescendant,
48)
49from camcops_server.cc_modules.cc_html import (
50 get_yes_no_none,
51 tr_qa,
52)
53from camcops_server.cc_modules.cc_request import CamcopsRequest
54from camcops_server.cc_modules.cc_sqla_coltypes import (
55 CamcopsColumn,
56 PendulumDateTimeAsIsoTextColType,
57)
58from camcops_server.cc_modules.cc_sqlalchemy import Base
59from camcops_server.cc_modules.cc_task import Task, TaskHasPatientMixin
60from camcops_server.cc_modules.cc_text import SS
62log = logging.getLogger(__name__)
65LOWER_MARKER = 0.25
66UPPER_MARKER = 0.75
67EQUATION_COMMENT = (
68 "logits: L(X) = intercept + slope * X; "
69 "probability: P = 1 / (1 + exp(-intercept - slope * X))"
70)
71MODALITY_AUDITORY = 0
72MODALITY_VISUAL = 1
73DP = 3
76# =============================================================================
77# CardinalExpDetThreshold
78# =============================================================================
80class CardinalExpDetThresholdTrial(GenericTabletRecordMixin, TaskDescendant,
81 Base):
82 __tablename__ = "cardinal_expdetthreshold_trials"
84 cardinal_expdetthreshold_id = Column(
85 "cardinal_expdetthreshold_id", Integer,
86 nullable=False,
87 comment="FK to CardinalExpDetThreshold"
88 )
89 trial = Column(
90 "trial", Integer,
91 nullable=False,
92 comment="Trial number (0-based)"
93 )
95 # Results
96 trial_ignoring_catch_trials = Column(
97 "trial_ignoring_catch_trials", Integer,
98 comment="Trial number, ignoring catch trials (0-based)"
99 )
100 target_presented = Column(
101 "target_presented", Integer,
102 comment="Target presented? (0 no, 1 yes)"
103 )
104 target_time = Column(
105 "target_time", PendulumDateTimeAsIsoTextColType,
106 comment="Target presentation time (ISO-8601)"
107 )
108 intensity = Column(
109 "intensity", Float,
110 comment="Target intensity (0.0-1.0)"
111 )
112 choice_time = Column(
113 "choice_time", PendulumDateTimeAsIsoTextColType,
114 comment="Time choice offered (ISO-8601)"
115 )
116 responded = Column(
117 "responded", Integer,
118 comment="Responded? (0 no, 1 yes)"
119 )
120 response_time = Column(
121 "response_time", PendulumDateTimeAsIsoTextColType,
122 comment="Time of response (ISO-8601)"
123 )
124 response_latency_ms = Column(
125 "response_latency_ms", Integer,
126 comment="Response latency (ms)"
127 )
128 yes = Column(
129 "yes", Integer,
130 comment="Subject chose YES? (0 didn't, 1 did)"
131 )
132 no = Column(
133 "no", Integer,
134 comment="Subject chose NO? (0 didn't, 1 did)"
135 )
136 caught_out_reset = Column(
137 "caught_out_reset", Integer,
138 comment="Caught out on catch trial, thus reset? (0 no, 1 yes)"
139 )
140 trial_num_in_calculation_sequence = Column(
141 "trial_num_in_calculation_sequence", Integer,
142 comment="Trial number as used for threshold calculation"
143 )
145 @classmethod
146 def get_html_table_header(cls) -> str:
147 return f"""
148 <table class="{CssClass.EXTRADETAIL}">
149 <tr>
150 <th>Trial# (0-based)</th>
151 <th>Trial# (ignoring catch trials) (0-based)</th>
152 <th>Target presented?</th>
153 <th>Target time</th>
154 <th>Intensity</th>
155 <th>Choice time</th>
156 <th>Responded?</th>
157 <th>Response time</th>
158 <th>Response latency (ms)</th>
159 <th>Yes?</th>
160 <th>No?</th>
161 <th>Caught out (and reset)?</th>
162 <th>Trial# in calculation sequence</th>
163 </tr>
164 """
166 def get_html_table_row(self) -> str:
167 return ("<tr>" + "<td>{}</td>" * 13 + "</th>").format(
168 self.trial,
169 self.trial_ignoring_catch_trials,
170 self.target_presented,
171 self.target_time,
172 ws.number_to_dp(self.intensity, DP),
173 self.choice_time,
174 self.responded,
175 self.response_time,
176 self.response_latency_ms,
177 self.yes,
178 self.no,
179 ws.webify(self.caught_out_reset),
180 ws.webify(self.trial_num_in_calculation_sequence)
181 )
183 # -------------------------------------------------------------------------
184 # TaskDescendant overrides
185 # -------------------------------------------------------------------------
187 @classmethod
188 def task_ancestor_class(cls) -> Optional[Type["Task"]]:
189 return CardinalExpDetThreshold
191 def task_ancestor(self) -> Optional["CardinalExpDetThreshold"]:
192 return CardinalExpDetThreshold.get_linked(
193 self.cardinal_expdetthreshold_id, self)
196class CardinalExpDetThreshold(TaskHasPatientMixin, Task):
197 """
198 Server implementation of the Cardinal_ExpDetThreshold task.
199 """
200 __tablename__ = "cardinal_expdetthreshold"
201 shortname = "Cardinal_ExpDetThreshold"
202 use_landscape_for_pdf = True
204 # Config
205 modality = Column(
206 "modality", Integer,
207 comment="Modality (0 auditory, 1 visual)"
208 )
209 target_number = Column(
210 "target_number", Integer,
211 comment="Target number (within available targets of that modality)"
212 )
213 background_filename = CamcopsColumn(
214 "background_filename", Text,
215 exempt_from_anonymisation=True,
216 comment="Filename of media used for background"
217 )
218 target_filename = CamcopsColumn(
219 "target_filename", Text,
220 exempt_from_anonymisation=True,
221 comment="Filename of media used for target"
222 )
223 visual_target_duration_s = Column(
224 "visual_target_duration_s", Float,
225 comment="Visual target duration (s)"
226 )
227 background_intensity = Column(
228 "background_intensity", Float,
229 comment="Intensity of background (0.0-1.0)"
230 )
231 start_intensity_min = Column(
232 "start_intensity_min", Float,
233 comment="Minimum starting intensity (0.0-1.0)"
234 )
235 start_intensity_max = Column(
236 "start_intensity_max", Float,
237 comment="Maximum starting intensity (0.0-1.0)"
238 )
239 initial_large_intensity_step = Column(
240 "initial_large_intensity_step", Float,
241 comment="Initial, large, intensity step (0.0-1.0)"
242 )
243 main_small_intensity_step = Column(
244 "main_small_intensity_step", Float,
245 comment="Main, small, intensity step (0.0-1.0)"
246 )
247 num_trials_in_main_sequence = Column(
248 "num_trials_in_main_sequence", Integer,
249 comment="Number of trials required in main sequence"
250 )
251 p_catch_trial = Column(
252 "p_catch_trial", Float,
253 comment="Probability of catch trial"
254 )
255 prompt = CamcopsColumn(
256 "prompt", UnicodeText,
257 exempt_from_anonymisation=True,
258 comment="Prompt given to subject"
259 )
260 iti_s = Column(
261 "iti_s", Float,
262 comment="Intertrial interval (s)"
263 )
265 # Results
266 finished = Column(
267 "finished", Integer,
268 comment="Subject finished successfully (0 no, 1 yes)"
269 )
270 intercept = Column(
271 "intercept", Float,
272 comment=EQUATION_COMMENT
273 )
274 slope = Column(
275 "slope", Float,
276 comment=EQUATION_COMMENT
277 )
278 k = Column(
279 "k", Float,
280 comment=EQUATION_COMMENT + "; k = slope"
281 )
282 theta = Column(
283 "theta", Float,
284 comment=EQUATION_COMMENT + "; theta = -intercept/k = -intercept/slope"
285 )
287 # Relationships
288 trials = ancillary_relationship(
289 parent_class_name="CardinalExpDetThreshold",
290 ancillary_class_name="CardinalExpDetThresholdTrial",
291 ancillary_fk_to_parent_attr_name="cardinal_expdetthreshold_id",
292 ancillary_order_by_attr_name="trial"
293 ) # type: List[CardinalExpDetThresholdTrial]
295 @staticmethod
296 def longname(req: "CamcopsRequest") -> str:
297 _ = req.gettext
298 return _("Cardinal RN – Threshold determination for "
299 "Expectation–Detection task")
301 def is_complete(self) -> bool:
302 return bool(self.finished)
304 def get_trial_html(self, req: CamcopsRequest) -> str:
305 """
306 Note re plotting markers without lines:
308 .. code-block:: python
310 import matplotlib.pyplot as plt
312 fig, ax = plt.subplots()
313 ax.plot([1, 2], [1, 2], marker="+", color="r", linestyle="-")
314 ax.plot([1, 2], [2, 1], marker="o", color="b", linestyle="None")
315 fig.savefig("test.png")
316 # ... the "absent" line does NOT "cut" the red one.
318 Args:
319 req:
321 Returns:
323 """
324 trialarray = self.trials
325 html = CardinalExpDetThresholdTrial.get_html_table_header()
326 for t in trialarray:
327 html += t.get_html_table_row()
328 html += """</table>"""
330 # Don't add figures if we're incomplete
331 if not self.is_complete():
332 return html
334 # Add figures
336 figsize = (
337 PlotDefaults.FULLWIDTH_PLOT_WIDTH / 2,
338 PlotDefaults.FULLWIDTH_PLOT_WIDTH / 2
339 )
340 jitter_step = 0.02
341 dp_to_consider_same_for_jitter = 3
342 y_extra_space = 0.1
343 x_extra_space = 0.02
344 trialfig = req.create_figure(figsize=figsize)
345 trialax = trialfig.add_subplot(MatplotlibConstants.WHOLE_PANEL)
346 notcalc_detected_x = []
347 notcalc_detected_y = []
348 notcalc_missed_x = []
349 notcalc_missed_y = []
350 calc_detected_x = []
351 calc_detected_y = []
352 calc_missed_x = []
353 calc_missed_y = []
354 catch_detected_x = []
355 catch_detected_y = []
356 catch_missed_x = []
357 catch_missed_y = []
358 all_x = []
359 all_y = []
360 for t in trialarray:
361 x = t.trial
362 y = t.intensity
363 all_x.append(x)
364 all_y.append(y)
365 if t.trial_num_in_calculation_sequence is not None:
366 if t.yes:
367 calc_detected_x.append(x)
368 calc_detected_y.append(y)
369 else:
370 calc_missed_x.append(x)
371 calc_missed_y.append(y)
372 elif t.target_presented:
373 if t.yes:
374 notcalc_detected_x.append(x)
375 notcalc_detected_y.append(y)
376 else:
377 notcalc_missed_x.append(x)
378 notcalc_missed_y.append(y)
379 else: # catch trial
380 if t.yes:
381 catch_detected_x.append(x)
382 catch_detected_y.append(y)
383 else:
384 catch_missed_x.append(x)
385 catch_missed_y.append(y)
386 trialax.plot(
387 all_x, all_y,
388 marker=MatplotlibConstants.MARKER_NONE,
389 color=MatplotlibConstants.COLOUR_GREY_50,
390 linestyle=MatplotlibConstants.LINESTYLE_SOLID,
391 label=None
392 )
393 trialax.plot(
394 notcalc_missed_x, notcalc_missed_y,
395 marker=MatplotlibConstants.MARKER_CIRCLE,
396 color=MatplotlibConstants.COLOUR_BLACK,
397 linestyle=MatplotlibConstants.LINESTYLE_NONE,
398 label="miss"
399 )
400 trialax.plot(
401 notcalc_detected_x, notcalc_detected_y,
402 marker=MatplotlibConstants.MARKER_PLUS,
403 color=MatplotlibConstants.COLOUR_BLACK,
404 linestyle=MatplotlibConstants.LINESTYLE_NONE,
405 label="hit"
406 )
407 trialax.plot(
408 calc_missed_x, calc_missed_y,
409 marker=MatplotlibConstants.MARKER_CIRCLE,
410 color=MatplotlibConstants.COLOUR_RED,
411 linestyle=MatplotlibConstants.LINESTYLE_NONE,
412 label="miss, scored"
413 )
414 trialax.plot(
415 calc_detected_x, calc_detected_y,
416 marker=MatplotlibConstants.MARKER_PLUS,
417 color=MatplotlibConstants.COLOUR_BLUE,
418 linestyle=MatplotlibConstants.LINESTYLE_NONE,
419 label="hit, scored"
420 )
421 trialax.plot(
422 catch_missed_x, catch_missed_y,
423 marker=MatplotlibConstants.MARKER_CIRCLE,
424 color=MatplotlibConstants.COLOUR_GREEN,
425 linestyle=MatplotlibConstants.LINESTYLE_NONE,
426 label="CR"
427 )
428 trialax.plot(
429 catch_detected_x, catch_detected_y,
430 marker=MatplotlibConstants.MARKER_STAR,
431 color=MatplotlibConstants.COLOUR_GREEN,
432 linestyle=MatplotlibConstants.LINESTYLE_NONE,
433 label="FA"
434 )
435 leg = trialax.legend(
436 numpoints=1,
437 fancybox=True, # for set_alpha (below)
438 loc="best", # bbox_to_anchor=(0.75, 1.05)
439 labelspacing=0,
440 handletextpad=0,
441 prop=req.fontprops
442 )
443 leg.get_frame().set_alpha(0.5)
444 trialax.set_xlabel("Trial number (0-based)", fontdict=req.fontdict)
445 trialax.set_ylabel("Intensity", fontdict=req.fontdict)
446 trialax.set_ylim(0 - y_extra_space, 1 + y_extra_space)
447 trialax.set_xlim(-0.5, len(trialarray) - 0.5)
448 req.set_figure_font_sizes(trialax)
450 fitfig = None
451 if self.k is not None and self.theta is not None:
452 fitfig = req.create_figure(figsize=figsize)
453 fitax = fitfig.add_subplot(MatplotlibConstants.WHOLE_PANEL)
454 detected_x = []
455 detected_x_approx = []
456 detected_y = []
457 missed_x = []
458 missed_x_approx = []
459 missed_y = []
460 all_x = []
461 for t in trialarray:
462 if t.trial_num_in_calculation_sequence is not None:
463 all_x.append(t.intensity)
464 approx_x = f"{t.intensity:.{dp_to_consider_same_for_jitter}f}" # noqa
465 if t.yes:
466 detected_y.append(
467 1 -
468 detected_x_approx.count(approx_x) * jitter_step)
469 detected_x.append(t.intensity)
470 detected_x_approx.append(approx_x)
471 else:
472 missed_y.append(
473 0 + missed_x_approx.count(approx_x) * jitter_step)
474 missed_x.append(t.intensity)
475 missed_x_approx.append(approx_x)
476 fit_x = np.arange(0.0 - x_extra_space, 1.0 + x_extra_space, 0.001)
477 fit_y = logistic(fit_x, self.k, self.theta)
478 fitax.plot(
479 fit_x, fit_y,
480 color=MatplotlibConstants.COLOUR_GREEN,
481 linestyle=MatplotlibConstants.LINESTYLE_SOLID
482 )
483 fitax.plot(
484 missed_x, missed_y,
485 marker=MatplotlibConstants.MARKER_CIRCLE,
486 color=MatplotlibConstants.COLOUR_RED,
487 linestyle=MatplotlibConstants.LINESTYLE_NONE
488 )
489 fitax.plot(
490 detected_x, detected_y,
491 marker=MatplotlibConstants.MARKER_PLUS,
492 color=MatplotlibConstants.COLOUR_BLUE,
493 linestyle=MatplotlibConstants.LINESTYLE_NONE
494 )
495 fitax.set_ylim(0 - y_extra_space, 1 + y_extra_space)
496 fitax.set_xlim(np.amin(all_x) - x_extra_space,
497 np.amax(all_x) + x_extra_space)
498 marker_points = []
499 for y in (LOWER_MARKER, 0.5, UPPER_MARKER):
500 x = inv_logistic(y, self.k, self.theta)
501 marker_points.append((x, y))
502 for p in marker_points:
503 fitax.plot(
504 [p[0], p[0]], # x
505 [-1, p[1]], # y
506 color=MatplotlibConstants.COLOUR_GREY_50,
507 linestyle=MatplotlibConstants.LINESTYLE_DOTTED
508 )
509 fitax.plot(
510 [-1, p[0]], # x
511 [p[1], p[1]], # y
512 color=MatplotlibConstants.COLOUR_GREY_50,
513 linestyle=MatplotlibConstants.LINESTYLE_DOTTED
514 )
515 fitax.set_xlabel("Intensity", fontdict=req.fontdict)
516 fitax.set_ylabel("Detected? (0=no, 1=yes; jittered)",
517 fontdict=req.fontdict)
518 req.set_figure_font_sizes(fitax)
520 html += f"""
521 <table class="{CssClass.NOBORDER}">
522 <tr>
523 <td class="{CssClass.NOBORDERPHOTO}">
524 {req.get_html_from_pyplot_figure(trialfig)}
525 </td>
526 <td class="{CssClass.NOBORDERPHOTO}">
527 {req.get_html_from_pyplot_figure(fitfig)}
528 </td>
529 </tr>
530 </table>
531 """
533 return html
535 def logistic_x_from_p(self, p: Optional[float]) -> Optional[float]:
536 try:
537 return (math.log(p / (1 - p)) - self.intercept) / self.slope
538 except (TypeError, ValueError):
539 return None
541 def get_task_html(self, req: CamcopsRequest) -> str:
542 if self.modality == MODALITY_AUDITORY:
543 modality = req.sstring(SS.AUDITORY)
544 elif self.modality == MODALITY_VISUAL:
545 modality = req.sstring(SS.VISUAL)
546 else:
547 modality = None
548 h = f"""
549 <div class="{CssClass.SUMMARY}">
550 <table class="{CssClass.SUMMARY}">
551 {self.get_is_complete_tr(req)}
552 </table>
553 </div>
554 <div class="{CssClass.EXPLANATION}">
555 The ExpDet-Threshold task measures visual and auditory
556 thresholds for stimuli on a noisy background, using a
557 single-interval up/down method. It is intended as a prequel to
558 the Expectation–Detection task.
559 </div>
560 <table class="{CssClass.TASKCONFIG}">
561 <tr>
562 <th width="50%">Configuration variable</th>
563 <th width="50%">Value</th>
564 </tr>
565 """
566 h += tr_qa("Modality", modality)
567 h += tr_qa("Target number", self.target_number)
568 h += tr_qa("Background filename", ws.webify(self.background_filename))
569 h += tr_qa("Background intensity", self.background_intensity)
570 h += tr_qa("Target filename", ws.webify(self.target_filename))
571 h += tr_qa("(For visual targets) Target duration (s)",
572 self.visual_target_duration_s)
573 h += tr_qa("Start intensity (minimum)", self.start_intensity_min)
574 h += tr_qa("Start intensity (maximum)", self.start_intensity_max)
575 h += tr_qa("Initial (large) intensity step",
576 self.initial_large_intensity_step)
577 h += tr_qa("Main (small) intensity step",
578 self.main_small_intensity_step)
579 h += tr_qa("Number of trials in main sequence",
580 self.num_trials_in_main_sequence)
581 h += tr_qa("Probability of a catch trial", self.p_catch_trial)
582 h += tr_qa("Prompt", self.prompt)
583 h += tr_qa("Intertrial interval (ITI) (s)", self.iti_s)
584 h += f"""
585 </table>
586 <table class="{CssClass.TASKDETAIL}">
587 <tr><th width="50%">Measure</th><th width="50%">Value</th></tr>
588 """
589 h += tr_qa("Finished?", get_yes_no_none(req, self.finished))
590 h += tr_qa("Logistic intercept",
591 ws.number_to_dp(self.intercept,
592 DP))
593 h += tr_qa("Logistic slope",
594 ws.number_to_dp(self.slope, DP))
595 h += tr_qa("Logistic k (= slope)",
596 ws.number_to_dp(self.k, DP))
597 h += tr_qa("Logistic theta (= –intercept/slope)",
598 ws.number_to_dp(self.theta, DP))
599 h += tr_qa(f"Intensity for {100 * LOWER_MARKER}% detection",
600 ws.number_to_dp(self.logistic_x_from_p(LOWER_MARKER),
601 DP))
602 h += tr_qa("Intensity for 50% detection",
603 ws.number_to_dp(self.theta, DP))
604 h += tr_qa(f"Intensity for {100 * UPPER_MARKER}% detection",
605 ws.number_to_dp(self.logistic_x_from_p(UPPER_MARKER), DP))
606 h += """
607 </table>
608 """
609 h += self.get_trial_html(req)
610 return h