Copula Modelling of Agitation-Sedation (A-S) in ICU: Threshold Analysis of Nurses’ Scores of A-S and Automated Drug Infusions by Protocol

Abstract

Pain management is increasingly recognised as a formal medical subspecialty worldwide. Empirical distributions of the nurses’ ratings of a patient’s pain and/or agitation levels and the administered dose of sedative are often positively skewed, and if the joint distribution is non-elliptical, then high nurses’ ratings of a patient’s agitation levels may not correspond to the true occurrences of patient’s agitation-sedation (A-S). Copulas are used to capture such nonlinear dependence between skewed distributions and check for the presence of lower (LT) and/or upper tail (UT) dependence between the nurses’ A-S rating and the automated sedation dose, thus finding thresholds and regions of mismatch between the nurse’s scores and automated sedation dose, thereby suggesting a possible way forward for an improved alerting system for over- or under-sedation. We find for LT dependence nurses tend to underestimate the patient’s agitation in the moderate agitation zone. In the mild agitation zone, nurses tend to assign a rating, that is, on average, 0.30 to 0.45 points lower than expected for the patient’s given agitation severity. For UT dependence in the moderate agitation zone, nurses tend to either moderately or strongly underestimate patient’s agitation, but in periods of severe agitation, nurses tend to overestimate a patient’s agitation. Our approach lends credence to augmenting conventional RASS and SAS agitation measures with semi-automated systems and identifying thresholds and regions of deviance for alerting increased risk.

1. Introduction

Pain management is becoming increasingly recognised as a formal medical subspecialty worldwide. Pain is the most common reason that patients come to the emergency department. Emergency nurses have an indispensable role in the management of this pain [1]. Sedation in the intensive care unit (ICU) is challenging, as both over- and under-sedation are detrimental [2]. Optimal sedation and analgesic strategies, combined with delirium management, are difficult when caring for critically ill patients. For sedation monitoring, the most widely used tools are the Richmond Agitation and Sedation Scale (RASS) [3] and Sedation Agitation Scale (SAS) [4]. Assessments using RASS and SAS are undertaken intermittently and traditionally rely on patients’ behavioural response to stimulation, which perturbs rest and sleep [5, 6, 7, 8]. Various studies have suggested that a non-stimulating method for “continuous” sedation monitoring may be beneficial and allow for more frequent assessment.

Indeed earlier, Rudge, Chase, and Shaw [9] discussed target-controlled infusion (TCI) systems to deliver drugs to maintain target plasma concentrations, using a pharmacokinetic model, shown to be feasible when anaesthesia is given over short periods of reduced consciousness and well-known pharmacology is invoked. Infusion systems that regulate the infusion rate to maintain target agitation levels, to regulate the primary metric for long-term sedation, are one approach to improving care in the ICU. The data analysed in this chapter pertains to the scenario and data type studied earlier by Hudson, Rudge and colleagues [9, 10, 11, 12, 13, 14].

These authors have suggested that assessing the severity of agitation is a challenging clinical problem as variability related to drug metabolism for each individual is often subjective. A multitude of previous studies also suggest assessment accuracy of the sedation quality conducted by nurses tend to suffer from subjectivity and lead to sub-optimal sedation [9, 10, 11, 12, 13, 14]. For example, it has been recommended by some authors to use lighter than deeper levels of sedation. And that sedation should be reviewed and adjusted regularly [5, 6, 7, 8]. Agitation management methods frequently rely on subjective agitation assessment the carers then select an appropriate infusion rate based upon their evaluation of these scales, experience, and intuition [3, 4, 5]. This approach usually leads to largely continuous infusions which lack a bolus-focused approach, commonly resulting in over or under-sedation.

The work of [9, 10, 11, 12, 13, 14] aimed to enhance feedback protocols for medical decision support systems and eventually automated sedation administration. A minimal differential equation model to predict or simulate each patient’s agitation-sedation status over time was presented in [9] for our ICU patients and was shown to capture patient agitation-sedation (A-S) dynamics. The use of quantitative modelling to enhance understanding of the agitation-sedation (A-S) system and the provision of an A-S simulation platform is one of the key tools in this area of patient critical care. A more refined A-S model, which utilised regression with an Epanechnikov kernel was formulated by [9]. A Bayesian approach using densities and wavelet shrinkage methods was later suggested by [10] to assess a previously derived deterministic, parametric A-S model [11], thus successfully challenging the practice of sedating ICU patients using continuous infusions.

Wavelets approaches [10, 11] were shown to provide reliable diagnostics and visualisation tools to assess A-S models, giving alternative metrics of A-S control to assess the validity of the earlier A-S deterministic models (Table 3 in [10]). This suite of wavelet metrics based on the discrete wavelet transform (DWT) established the value of earlier deterministic agitation-sedation (A-S) models against empirical (recorded) dynamic A-S infusion profiles, providing robust performance metrics of A-S control and excellent tools, based on the classification of patients into poor and good trackers based on Wavelet Probability Bands (WPBs). Importantly, the WPBs were shown to be a useful patient-specific method by which to identify and detect regions in the patient’s A-S profile i.e., times whilst in ICU, where the simulated infusion rate performs poorly, thus providing visual and quantified ways to help improve and distil the deterministic A-S model and in practice be a gauge to alert carers.

The aim of this chapter is to identify regions of poor and good control using copulas. Copulas are functions that join or connect multivariate distribution functions to their one-dimensional marginal distribution functions. Copulas have had applications in fields such as finance [15, 16], public health and medicine [17], and actuarial science [18, 19]. Empirical distributions of the nurses’ ratings of a patient’s pain and/or agitation levels and the administered dose of sedative are often positively skewed and if the joint distribution is non-elliptical, then high nurses’ ratings of a patient’s agitation levels may not correspond to the occurrences of patient’s A-S profile with large infusion dose. Copulas are used as they measure nonlinear dependencies capturing the dependence between skewed distributions. Copulas are widely applied in diverse fields, including health services research and medical studies, quantitative risk management, econometric modelling, environmental studies, finance, and hydrology.

Advantages of using copulas in modelling are: (i) capacity to model both linear and non-linear dependence; (ii) allowing an arbitrary choice of a marginal distribution; and (iii) capability of modelling extreme endpoints. Copulas are functions that “couple together” the marginal cumulative distribution functions (CDFs) of a random vector to form its joint CDF. When used in statistical modelling, copulas can estimate multivariate distributions of data involving two or more outcome variables for mixed type, complex data. We determine the best-fit copula type for all patients with a focus on differences between poor and good trackers, where classification of patients into poor and good trackers was based on Wavelet Probability Bands (WPB) [10, 11].

This chapter builds on the earlier pilot work of Tursunalieva et al. [20, 21]] to address the gap in the methodology by integrating non-elliptical dependence structure between nurses’ rating of a patient’s agitation level and the automated sedation dose. In an earlier pilot work discussed by Hudson [22], the tail thresholds of two (2) test patients were determined manually, whereas in [21] the dynamic programming algorithm of Bai and Perron [23] was used to establish the lower and upper tail threshold. Copula mathematics allows us to determine and identify lower and/or upper tail thresholds when they exist for all 36 intensive care unit patients’ agitation-sedation profiles collected at Christchurch Hospital, School of Medicine and Health Sciences, NZ and analysed earlier in [9, 10, 11, 12, 13, 14]. Infusion data were recorded using an electronic drug infusion device for all admitted ICU patients during a nine-month observation period and required more than 24 hours of sedation.

In this chapter, our novel and general formulation of the equation relating each patient’s nurses’ score to the automated infusion dose is given by the following expression nurses’ score = intercept+α*Dose–β*Dose*LT region+ γ*Dose*UT region. This formulation accounts for the non-linear relationships between the nurses’ A-S rating and the automated sedation dose, and permits identification of thresholds and regions of mismatch between the nurse’s scores and sedation dose, thereby suggesting a possible way forward for an improved alerting system for over/under-sedation.

Establishing the presence of tail dependence and patient-specific thresholds for areas with different agitation intensities has significant implications for the effective administration of sedatives. Better management of A-S states will allow clinicians to improve the efficacy of care and reduce healthcare costs. Our approach lends credence to augmenting conventional RASS and SAS agitation measures with semi-automated systems [24, 25, 26] and identifying thresholds and regions of deviance for alerting increased risk. Better management of A-S states will allow clinicians to improve the efficacy of care and reduce healthcare costs.

1.1 Data and methods

This chapter models the agitation-sedation profiles of 36 patients collected at the Christchurch Hospital, Christchurch School of Medicine and Health Sciences, NZ. Two measures were recorded for each patient: (1) the nurses’ ratings/scores of a patient’s agitation level, and (2) an automated sedation dose (see Figure 1). Infusion data were recorded using an electronic drug infusion device for all admitted ICU patients during a nine-month observation period and required more than 24 hours of sedation. Infusion data containing less than 48 hours of continuous data, or data from patients whose sedation requirements were extreme, such as those with severe head injuries, were excluded [9, 10].

Figure 1.

Diagram of the feedback loop employing nursing staff’s feedback of subjectively assessed patient agitation through the infusion controller (diagram is sourced from Chase et al.) [12, 14].

A total of 36 ICU patients met these requirements and were enrolled in the study. Classification of patients into poor and good trackers, based on the Wavelet Probability Bands (WPB), is given in Table 1. The so-called good tracker delineates the scenario where the nurse’s rating scores remains within the (time-based) 90% coverage of wavelet probability band (WPB) based on the simulated dose profiles [10, 11]. Poor tracking delineates the scenario where the nurses’ rating scores remain outside the (time based) 90% coverage of wavelet probability band (WPB) for a significant portion of time based on the simulated dose profiles [11].

WPB [11]	WCORR [10]	Chase et al. [14]	Rudge et al. [9]
2	2	—	—
4	4	—	—
—	—	6	—
7	7	7	7
10	10	—	10
11	11	—	11
—	—	12	—
—	—	—	13
—	—	17	—
21	21	21	—
22	22	—	22
27	27	27	27
28	28	—	28
—	29	—	29
32	32	—	—
33	33	—	33
34	34	34	—
—	35	—	35
Total: N1 = 12	Total: N2 = 14	Total: N3 = 7	Total: N4 = 10

Table 1.

Patient numbers of the poor trackers according to the criteria of 4 studies, developed earlier in [11].

Low WPB 90% values under 70% indicate a poor tracker by Kang’s WPB diagnostics [11].

By way of illustration, we carefully examine four patients from the pool of 36 patients. Tables 2 and 3 summarise each of these 4 patients’ WPB tracker status, time to first, second and third violation outside the WPB bands, their total number of violations over ICU stay, and patient’s time in ICU, along with their specific WPB% value. Display of their line profiles of nurses’ rating of A-S in relation to drug infusion dose over time, for each of the 4 patients (P8, P27, P18, P28) are given in Figures 2–4. Note that a violation event occurs when the nurses observed agitation score or rating is outside either the lower or upper limits of the 90%WPB bands associated with the patient’s automated infusion dose trajectory over time in ICU.

	V1	V2	V3	Total V’s	Time in ICU	WPB%
P8/Good	1	2	3	46	128	87.5%
P27/Poor	1	4	5	89	225	43.7%

Table 2.

Time to the patient-specific, 1st violation denoted by V1, second violation V2 and third violation V3, total number of violations, total ICU time and WPB% values.

	V1	V2	V3	Total V’s	Time in ICU	WPB%
P18/Good	2	24	26	20	64	93.8%
P28/Poor	1	5	12	114	203	50.8%

Table 3.

Time to the patient-specific, 1st violation V1, second violation V2, and third violation V3, total number of violations, total ICU time and WPB% values.

Figure 2.

Line plot of nurses’ rating of patient agitation and the automated sedation dose for patient 8 (WPB-based good tracker).

Figure 3.

Line plot of nurses’ rating of patient agitation and the automated sedation dose for patient 27 (WPB-based poor tracker).

Figure 4.

Line plot (WPB% bands), patient 18 (LHS, good) and 28 (RHS, poor tracker).

The first patient (patient 8) in Table 2 is a good WPB tracker and the second (patient 27) a poor WPB tracker both studied in [20, 21], for which the upper tail thresholds of the nurses’ scores using copulas were established, in a pilot study of these 2 patients using copula mathematics. We also refer the reader to Hudson’s chapter in this book “Modelling Agitation-Sedation (A-S) in ICU: an Empirical Transition and Time to event analysis of poor and good tracking between nurses scores and automated A-S measures” [22].

The corresponding WPB% values for patient 8 and patient 27 are 87.5% and 43.7%, respectively (Table 2). Overall, the minimum, median and maximum WPB% values for the 24 good trackers is (58.8%, 87.5%, 96.9%) and (47.3%, 64.8%, 77.3%) for the 12 poor trackers (Table 1). Noteworthy also is that the A-S time series of these two patients examined (P8 and P27) were of disparate lengths - patient 8 had 10,561 time points and patient 27, had 13,441 time points. The full set of patients studied had a range of [3001-25,261] time points of automated dose assessments.

Patient 18 (good tracker) with a WPB% of 93.8% and patient 28 (poor tracker) with WPB% of 50.8% (Table 3) were studied in detail in [20, 21], for which both upper and lower tails/thresholds of over or under-estimation of agitation levels by the nurses’ rating were established using copula dependence analytics (see also [27, 28, 29]).

Patients vary according to their length of stay in ICU and consequently differ in their opportunity for violations to occur. The good trackers generally have shorter ICU time and thus less chance to exhibit an increased total number of violations.

The total number of WPB-based violations is greater for the poor trackers than for the good trackers, and it is the poor trackers that tend to have longer ICU times. There are three approximate categories of patient ICU time: 50−64, 113−128, and 205–256, and 19 of the 36 patients have an ICU time of ≤64.

2. Methodology

2.1 Copula formulation

In our study, the copula models aim to capture the dependence between the observed/recorded nurses’ rating and the automated sedation dose for each patient. We test for and utilise the so-called best fitting copula found. This section acts as a guide to decide if there exists a tail relationship between the nurses’ rating of patient agitation A-S score and automated sedation dose.

Analytically and contextually, the lower tail region corresponds to the patient specific mild agitation range and the upper tail region corresponds to the severe agitation range, with the non-tail middle region capturing the patient’s moderate agitation range. Clearly patient’s transition between these so-called mild, moderate and severe ranges or states over time in ICU. In our context, if the two distributions that is of nurses’ (observed) rating and that of the automated sedation dose are independent univariate Gaussians, then we can define the multivariate Gaussian distribution as the best fit. Let X and Y be independent Gaussian (with arbitrary means and variances), then ${\mathit{Z}}_{\mathit{j}}={\mathit{a}}_{\mathit{j}\mathbf{1}}\mathit{X}+{\mathit{a}}_{\mathit{j}\mathbf{2}}\mathit{Y}$ is univariate Gaussian for j = 1, 2, …, n and ${\mathit{a}}_{\mathit{j}\mathbf{1}},{\mathit{a}}_{\mathit{j}\mathbf{2}}$ are real constants, and Z is multivariate Gaussian.

Let$\mathbf{Z}=\left(\begin{array}{c}{\mathit{Z}}_{\mathbf{1}}\\ ⋮\\ {\mathit{Z}}_{\mathit{n}}\end{array}\right)=\mathit{A}\left(\begin{array}{c}\mathit{X}\\ \mathit{Y}\end{array}\right)$ be constructed based on the $n\times n$matrix $A$ where X and Y.

are independent and identically distributed (i.i.d.) N(0,1) random variables. Then $\mathrm{Z}+\mathrm{\mu }$ is multivariate Gaussian with mean vector $\mathrm{\mu }$ and covariance matrix $=A{A}^T$. From the Central Limit Theorem, the Gaussian distribution arises as a limit of a scaled sum of weakly dependent random variables [27, 28].

Parametric copula families are conventionally constructed to satisfy different combinations of bivariate dependence structures with tail behaviours [28]. The general definition of Archimedean copulas is given by Boateng [29]. It is noteworthy that the Clayton, Gumbel, and Frank copulas are examples of existing Archimedean copulas. A discussion about the Clayton, Gumbel, and Frank copulas, and tail dependence of a bivariate copula, Kendall’s tau representations, and copula models of the Clayton, Frank, and Gumbel copulas are also defined in [29]. The Clayton copula, for example, accommodates only lower tail dependence [30], the Frank copula allows dependence around the mode [31], and the Gumbel is relevant only when upper tail dependence exists [32]. The difference between the Clayton and Gumbel copulas is: (i) for the Clayton copula, correlations on the extreme left sides of distributions are more concentrated (i.e., higher correlations) than those in the extreme right sides of the distributions, and (ii) for the Gumbel copula, the correlations on the extreme right sides of distributions are more concentrated (i.e., higher correlations) than those in the extreme left sides of the distributions. The visuals in Section 3 of this chapter illustrate these trends. We refer the reader to Boeting [29] and below give the bivariate Gaussian formulation and bivariate Frank and Gumbel copula, as three examples.

Bivariate Gaussian copula

The copula cdf is:

$$\mathit{C}\left(\mathit{u},\mathit{\upsilon },\mathit{\rho }\right)={\mathbf{\Phi }}_{\mathbf{2}}\left({\mathbf{\Phi }}^{-\mathbf{1}}\left(\mathit{u}\right),{\mathbf{\Phi }}^{-\mathbf{1}}\left(\mathit{\upsilon }\right);\mathit{\rho }\right),\mathbf{0}<\mathit{u},\mathit{\upsilon }<\mathbf{1}.$$E1

Bivariate Frank

For $-\infty <\delta <\infty$, the copula cdf is:

$$\mathit{C}\left(\mathit{u},\mathit{\upsilon },\mathit{\delta }\right)=-{\mathit{\delta }}^{-\mathbf{1}}\mathit{log}\left(\frac{\mathbf{1}-{\mathit{e}}^{-\mathit{\delta }}-\left(\mathbf{1}-{\mathit{e}}^{-\mathit{\delta u}}\right)\left(\mathbf{1}-{\mathit{e}}^{-\mathit{\delta \upsilon }}\right)}{\mathbf{1}-{\mathit{e}}^{-\mathit{\delta }}}\right),\mathbf{0}<\mathit{u},\mathit{\upsilon }<\mathbf{1}.$$E2

Bivariate Gumbel copula

The copula cdf is: $\mathit{C}\left(\mathit{u},\mathit{\upsilon },\mathit{\delta }\right)=\mathit{exp}\left\{-{\left({\left[-\mathbf{log}\mathit{u}\right]}^{\mathit{\delta }}+{\left[-\mathbf{log}\mathit{\upsilon }\right]}^{\mathit{\delta }}\right)}^{\mathbf{1}/\mathit{\delta }}\right\}$,$\mathbf{0}\le \mathit{u},\mathit{\upsilon }\le \mathbf{1},\mathbf{1}\le \mathit{\delta }\le \infty .$ Upper tail dependence function for Gumbel copula is:

$${\mathit{b}}_{\mathit{U}}\left({\mathit{w}}_{\mathbf{1}}{\mathit{w}}_{\mathbf{2}};\mathit{\delta }\right)={\mathit{w}}_{\mathbf{1}}+{\mathit{w}}_{\mathbf{2}}-{\left({\mathit{w}}_{\mathbf{1}}^{\mathit{\delta }}+{\mathit{w}}_{\mathbf{2}}^{\mathit{\delta }}\right)}^{\mathbf{1}/\mathit{\delta }}.$$E3

2.2 Kendall K-plot construction

The best fitting copula was selected by maximum likelihood estimation, except for the t-copula, for which the degrees of freedom parameter is found by a crude profile likelihood optimisation over the interval (2, 10]. We use the Kendall plot (K-plot) [33] to determine the bivariate patient-specific thresholds in the cases where the best fitting copula has tails. The K-plot splits the data into two regions with significantly different strengths of dependence between nurses’ rating and the automated sedation dose, namely: (1) the main region with an approximately linear relationship; and (2) the tail regions with a non-linear relationship. Recently the K-plot has gained popularity with regard to its association with the receiver operating characteristic (ROC) curve, a pivotal biostatistical graphical tool traditionally used for testing the ability of biomarkers to discriminate between populations [34].

The K-plot adopts the familiar probability plot (Q-Q plot) to detect dependence. A lack of linearity of the standard Q-Q plot is an indication of non-normality of the distribution of a random variable. Similarly, in the absence of association between two variables, the K-plot is close to a straight line, while the amount of curvature in the K-plot is characteristic of the degree of dependence in the data, and is related, in a definite way, to the underlying copula. This method is closely related to Kendall’s tau statistic [35] from which it takes the name. For more details refer to [27, 33].

To construct a K-plot, we need to compute $H_i$ defined for a given pair ($X_i,Y_i)$with $1\le i\le n$ as follows: $H_i=\#\left\{j\ne i:X_j\le X_i,Y_j\le Y_i\right\}/\left(n-1\right)$. Next, we need to order the variable $H_i$, $H_{\left(1\right)}\le \dots \le H_{\left(n\right)}$ and plot the pairs $\left(W_{1:n},H_{\left(i\right)}\right),1\le i\le n,$where W_(1,n) is the expectation of the ith order statistic in a random sample of size n from the distribution $K_0$ of the $H_i$ under the null hypothesis of independence. Using the definition of the density of an order statistic, we define the form of $K_0$ under the null hypothesis of the independence, as follows:

$$W_{1:n}=n\left(\begin{array}{c}n-1\\ i-1\end{array}\right){\int }_0^1\omega {\left\{K_0\left(\omega \right)\right\}}^{i-1}\times {\left\{1-K_0\left(\omega \right)\right\}}^{n-1}d{K}_0\left(\omega \right),1\le i\le n.$$E4

2.3 Multiple threshold identification via dynamic programming algorithm

To identify a patient-specific threshold, we apply the dynamic programming algorithm discussed in Section 3.3 of [23] to use the dependence measure $H_i$. This algorithm captures multiple thresholds; however, to be consistent with the objective of this paper, we focus on determining the lower and upper tail thresholds. In our study, the lower tail threshold corresponds to the lowest (lower) threshold and the upper tail threshold corresponds to the highest (upper) threshold, respectively for either dose or the nurses’ score profiles.

2.4 Prediction equation of nurses’ score with respect to dose allowing for tails

Below we detail, as an illustration to the novel method that accommodates lower and upper thresholds beyond conventional correlational analysis using Kendall tau and copulas, the resultant equations specific to two patients, Patient 20 and Patient 8, which are both good trackers. P20 has no tails and P8 a lower tail. All the patients’ equations and their details are tabulated in the Appendix A.

The general formulation of the equation relating each patient’s nurses’ score to the automated infusion dose is given by either of the following expressions:

• $\begin{array}{l}\text{Score}=\text{intercept}+{\mathrm{\alpha }}^{\ast }\text{Dose}–\mathrm{\beta }{\text{Dose}}^{\ast }\mathrm{LT}+\mathrm{\gamma }{\text{Dose}}^{\ast }\mathrm{UT}\end{array}$

• $\begin{array}{l}\text{Score}=\mathrm{i}\text{ntercept}+\text{slope(Dose)}–\left(\text{slope}:\text{Lower region}\right)\text{Dose}\\ +\left(\text{slope dose}:\text{Upper region}\right).\end{array}$

Patient 20: The recorded or so-called nurses’ observed A-S score and the automated sedation dose for patient 20 are independent univariate Gaussians, therefore, their joint distribution is Gaussian (Tables 4 and 5 and LHS of Figure 5). This gives a bivariate Gaussian copula, with neither lower nor upper tails (RHS of Figure 5). Thus, the relationship between P20’s nurses’ score and infusion dose is estimated using a simple linear regression (SLR) equation, score = intercept + α*dose, with intercept and slope parameters − 0.31 and α = 1.16, respectively (Table 6 and RHS of Figure 5).

Patient 20 – good, no tails	Estimate	Standard Error	t value	P value
Intercept	−0.311	0.22	−1.41	0.16	NS
Dose	1.156	0.104	11.06	<2e-16	***

Table 4.

Patient 20 equation components and p-values, good tracker, no tails.

Adjusted R-squared: 0.49, Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘’ 1.

Patient 8 - good, lower tail	Estimate	Standard Error	t value	P value
Dose	1.012	0.0848	14.3	<2e-16	***
Dose*LT	−0.2849	0.0966	−2.95	0.0045	**

Table 5.

Patient 8 (P8) equation components - good tracker, lower tail tau = 0.70.

Adjusted R-squared: 0.766, Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘’ 1.

Figure 5.

Bivariate plots (LHS), copula type, and relevant tau (RHS) relating P 20’s nurses’ A-S score with dose, good tracker, no tails. SLR line is not shown.

Copula types	Tail	No. of copulas	WPB good trackers (Patient no.)	WPB poor trackers (Patient no.)
Clayton	Lower	2	P3, P19
Rotated Tawn type 1, 180 degrees	Lower	1		P35
Survival Gumbel	Lower	7	P1, P5, P6, P8, P14	P2, P7
Survival Joe	Lower	2	P15	P28
BB8	No tails	2	P30, P36
Frank	No tails	5	P13	P10, P11, P27, P33
Gaussian	No tails	5	P12, P17, P20, P29	P21
Survival BB8	No tails	7	P16, P18, P24, P26, P31	P22, P34
Survival BB7	Two tails	1		P32
t copula	Two Tails	1	P25
Gumbel	Upper	2	P23	P4
Joe	Upper	1	P37

Table 6.

List of copulas selected as optimal for the 36 patients (13 poor trackers). Shaded rows indicate the poor trackers. PX denotes patient number X.

Patient 8: In contrast to patient 20, for patient 8 (good tracker) the nurses’ recorded score and the automated sedation dose are skewed distributions (Figure 6), with the joint distribution being Survival Gumbel with a lower tail dependence. (lower tail tau = 0.70) (see Table 7 and Figure 6, top panel). Hence, the relationship between P8’s nurses’ score and the dose is estimated by our novel prediction equation, as follows.

Figure 6.

P8 bivariate plots (LHS), copula type and main region tau (upper panel) and K-plot and lower tail (LT) tau = 0.70, with fit lines (bottom panel) relating P8’s score with dose. P8 is a good tracker with lower tail thresholds for dose and nurses’ score: LT dose threshold = 3.02 and LT nurses’ score threshold = 2.57.

Patient no.	WPB status	Time			Copula tau (τ)
		n	Copula type	Tail	Main region	Lower tail	Upper tail
1		59	Survival Gumbel	Lower	0.56	0.67
2	Poor	63	Survival Gumbel	Lower	0.53	0.64
3		63	Clayton	Lower	0.49	0.66
4	Poor	63	Gumbel	Upper	0.61		0.69
5		48	Survival Gumbel	Lower	0.44	0.55
6		63	Survival Gumbel	Lower	0.48	0.57
7	Poor	63	Survival Gumbel	Lower	0.82	0.85
8		127	Survival Gumbel	Lower	0.63	0.70
10	Poor	255	Frank		0.67
11	Poor	111	Frank		0.71
12		127	Gaussian		0.56
13		63	Frank		0.74
14		49	Survival Gumbel	Lower	0.65	0.72
15		63	Survival Joe	Lower	0.53	0.74
16		127	Survival BB8		0.67
17		63	Gaussian		0.52
18		63	Survival BB8		0.69
19		127	Clayton	Lower	0.54	0.69
20		127	Gaussian		0.58
21	Poor	61	Gaussian		0.57
22	Poor	127	Survival BB8		0.64
23		57	Gumbel	Upper	0.61		0.68
24		127	Survival BB8		0.64
25		63	t copula	Both	0.57	0.06	0.06
26		63	Survival BB8		0.64
27	Poor	223	Frank		0.70
28	Poor	203	Survival Joe	Lower	0.59	0.79
29		53	Gaussian		0.68
30		60	BB8		0.46
31		255	Survival BB8		0.68
32	Poor	252	Survival BB7	Both	0.58	0.73	0.40
33	Poor	255	Frank		0.67
34	Poor	127	Survival BB8		0.60
35	Poor	211	Rot Tawn 1180 ^0	Lower	0.67	0.74
36		57	BB8		0.61
37		123	Joe	Upper	0.56		0.77

Table 7.

List of optimal copulas for the 36 patients with each patient’s WPB status, number of time points over ICU stay, copula type, values for copula tau (τ) for the main region, the lower or upper regions taus for the patient’s (score, dose) bivariate profiles.

Shaded rows indicate the poor trackers.

Score = 1.01Dose-0.29Dose*LT: The corresponding slope and lower tail parameters are α = 1.01 and β = −0.29 (Table 7). The highly significant slope of 1.01 (p < 0.00001) indicates that in this patient’s moderate agitation zone nurses tend to estimate agitation severity quite accurately. However, in the mild agitation zone, nurses tend to assign a rating that is, on average, 0.29 points lower than expected for the patient’s given (by automated dose) agitation severity (p = 0.0045) (RHS of bottom panel of Figure 6). There are 33 (out of a total of 127) such occurrences indicating that in approximately one in every four ratings, nurses tend to underestimate this patients’ agitation severity, LT dose threshold = 3.02 and LT score threshold = 2.57.

The Kendall-plot is then used to identify the lower tail threshold which occurs at the 26th percentile in the patient’s bivariate (dose, score) trajectory (Table A.1 Appendix), estimated by the algorithm of Bai and Perron [23]. For patient 8, the LT infusion dose threshold is 3.02 and LT nurse’s score threshold is 2.57 (see percentiles in Table A.1) and the full set of equations in Table A.2 in Appendix.

3. Results

3.1 Copula types across WPB status and tails status, and copula taus

Table 4 Gives the summary of copulas distributions selected as optimal for the total 36 patients stratified by WPB tracking status, copula type, number and type of tail(s), as gleaned from the more detailed Table 5, which gives each patient’s WPB status, number of time points captured over ICU stay, patient’s copula type, in addition to the values for copula tau (τ) for the main region, the lower or upper regions, as applicable for the patients’ (nurses’ score, dose) bivariate profiles.

In regards to the distribution of copula types by WPB (good versus poor) tracking status, there are four Frank copula types of the 13 poor trackers, in comparison to 1 good tracker being a Frank copula type. Of the 23 good trackers, the majority are either Survival Gumbel (5 of these), Survival BB8 (5 such), or Gaussian (4) copula types (Table 4).

Furthermore, of the 23 good trackers, eight patients (35%) have bivariate dependence with a lower tail, 12 patients (52%) with no tails, 1 patient P32 (4.3%) with both upper and lower tails, and 2 (8.7%) (P23 and 27) with an upper tail. In comparison of the 13 poor trackers, there are four patients (31%) with a lower tail, seven (54%) with no tails, one patient P25 (7.7%) with both tails, and one patient 4 (7.7%) with an upper tail (Table 4).

Copulas that are unique to the poor trackers are the Survival BB7 (P32 with upper and lower tails, a poor tracker, displayed in Figure 7 and the Rotated Tawn type 1, 180⁰ (P35 with an upper tail, a poor tracker displayed in Figure 8). Unique copula types found only for the good trackers are the Clayton (P3 and P19, each have lower tails), the BB8 (P30 and P36 both with no tails), the t copula (P25 with upper and lower tails) shown in Figure 9 and the Joe copula (P37, upper tail) displayed in Figure 10.

Figure 7.

P32 a poor tracker with 2 tails - copula plot, main tau, K- plot and tail taus (upper panel); best fit line(s) (lower panel). LT tau = 0.73, UT tau = 0.40, (LT, UT) dose threshold = (3.90, 7.02), (LT, UT) score threshold = (3.62, 7.63).

Figure 8.

P28 poor tracker with a lower tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P28 nurses’ score with dose. LT tail tau = 0.79, LT dose threshold = 1.70 and LT score threshold = 1.20.

Figure 9.

P25 good tracker with 2 tails - copula plot, main region tau K-plot, and tail taus (upper panel); best fit line(s) (lower panel) LT tau = 0.06, UT tau = 0.06, and (LT, UT) dose threshold = (2.68, 4.43), with (LT, UT) score threshold = (2.41, 3.94).

Figure 10.

P37 good tracker with upper tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P37 nurses’ A-S score with dose. UT tau = 0.77, UT dose threshold = 2.53, and UT nurse score threshold = 2.99.

3.2 Visuals, tail dependence, taus, dose/score thresholds of 10 patients

Details of equations, copulas and visualisations will focus on the following list of 10 patients given in Tables 8 and 9. Visualisations comprise patients’ score and dose trajectories with associated 95% WPB bands, copula plots, K-plots, along with associated equations, relevant tau for tails and a clear delineation dose and nurses’ score thresholds related to upper and lower tails, when they exist. Specifically lower and upper tail tau (τ) values, their percentile positions in the patient bivariate trajectories of length in ICU, associated lower and/or upper tail thresholds for the infusion dose and nurse’s score threshold are given in Table 8. Table 9 gives the associated equation parameters contrasting linear correlation (r) with our novel copula-tail based approach. Figures 11 and 12 display the bivariate trajectories and 90% WPB bands of the patients. Full set of equations and tail thresholds are given in Tables A.1 and A.2 in the Appendix.

			Lower Tail tau (τ)		Upper Tail tau (τ)		Dose threshold		Score threshold
P no.	Tails	n	n	Per centile	n	Per centile	Lower	Upper	Lower	Upper
2	Lower	63	28	44.44			0.96		0.83
4	Upper	63			15	76.19		2.82		2.75
15	Lower	63	22	34.92			3.47		3.05
23	Upper	57			10	82.46		1.77		1.97
25	Both tails	63	14	33.33	15	68.25	2.68	4.43	2.41	3.94
27	No tails	223
28	Lower	203	59	29.06			1.70		1.20
32	Both tails	252	44	20.64	38	74.60	3.90	7.02	3.62	7.63
35	Lower	211	42	19.91			0.72		0.66
37	Upper	123			17	86.18		2.53		2.99

Table 8.

Upper and/or lower tail positions (percentiles) and associated lower/upper thresholds for dose and nurses’ A-S score.

Shaded rows are the poor trackers.

		Linear correlation (r)			Novel Regression (estimates)				adj R^2 quare
P no.	Tail status	Main region	Lower tail	Upper tail	Intercept	Dose Slope	Dose Slope Lower	Dose Slope Upper	Simple LR R^2	Novel method R^2
2	Lower	0.39	0.77		0.27	0.85	−0.37		0.42	0.44
4	Upper	0.38		0.66	0.79	0.52		0.45	0.57	0.61
15	Lower	0.28	0.75		0.86	0.83	−0.34		0.45	0.47
23	Upper	0.51		0.62	0.34	0.58		0.52	0.64	0.68
25	Both	0.44	0.41	0.00	1.55	0.46	−0.35	0.34	0.60	0.70
27	None	0.84			−0.36	1.06			0.71
28	Lower	0.46	0.84		0.008	0.97	−0.35		0.57	0.59
32	Both	0.26	0.55	0.61	3.23	0.43	−0.75	0.33	0.60	0.67
35	Lower	0.72	0.77		0.23	0.91	−0.79		0.63	0.64
37	Upper	0.54		0.75	0.17	0.78		0.32	0.74	0.75

Table 9.

Upper and/or lower tail dose relationships, equations, and associated change in adjusted R² for simple LR vs our novel approach. Shaded rows indicate the poor trackers.

Figure 11.

Line plots of nurses’ score (observed, red) and dose (black line), with 95% WPB bands for P32, P25 (both tails); P4, P23 (upper tail); P2, P15 (lower tail).

Figure 12.

Line plots of nurses’ score (observed, red) vs. dose (black line), with 90% WPB bands for P37 good tracker.

3.3 Prediction equation of score in relation to dose allowing tail dependence

In this subsection, we focus on 7 of the 10 selected patients reported in Table 8. Specifically, two patients with both lower and upper tails (P32 and P35), three patients with only upper tails (P4, P23, P37), and two patients with only lower tails (P2 and P15). For each of these seven patients, our novel equation relating each patient’s nurses’ agitation severity rank score versus the patient’s infusion dose with parameter estimates and p-values is reported below. In addition, interpretation of the equations in regard to regions where the nurses’ score either under or overestimates the patient’s agitation with respect to so-called ground truth, this being the patient’s automated infusion dose is reported per patient.

Patient 32 poor tracker, 2 tails, R² squared (non adj, adjusted) = (0.60, 0.67).

Score = 3.23 + 0.43Dose-0.75Dose*LT + 0.33Dose*UT (Table 10): The intercept of 3.23 indicates that the patient is experiencing severe “chronic” background agitation. The slope of 0.43 indicates that in the moderate agitation zone nurses tend to strongly underestimate the patient’s agitation severity. In mild agitation zone, nurses tend to still underestimate the agitation severity. In severe agitation zone, nurses tend to overestimate the patient’s agitation severity on average 0.33 points higher (Figures 7 and 11). For P32 LT tau = 0.73, UT tau = 0.40, and its bivariate (LT, UT) dose thresholds are (3.90, 7.02), with (LT, UT) nurses’ score thresholds of (3.62, 7.63).

Patient 32 - poor, 2 tails Adjusted R-squared: 0.66			Estimate		Standard Error	t value		P value
Intercept			3.2342		0.5162	6.27		1.60E−09	***
Dose			0.4285		0.0828	5.18		4.70E−07	***
Dose*LT			−0.7506		0.1503	−4.99		1.10E−06	***
Dose*UT			0.3309		0.0497	6.65		1.80E−10	***
Patient 25 - good, 2 tails; Adjusted R-squared: 0.70			Estimate		Standard Error	t value		P value
Intercept			1.5481		0.5503	2.81		0.0067	**
Dose			0.4603		0.1543	2.98		0.0042	**
Dose*LT			−0.3475		0.1719	−2.02		0.0477	*
Dose*UT			0.3424		0.0788	4.34		5.60E-05	***
		Linear correlation (r)			Novel Regression (estimates)				adj R^2quare
P no.	Tail status	Main region	Lower tail	Upper tail	Intercept	Dose Slope	Dose Slope Lower	Dose Slope Upper	Simple LR R^2	Novel method R^2
32	Both	0.26	0.55	0.61	3.23	0.43	−0.75	0.33	0.60	0.67
25	Both	0.44	0.41	0.00	1.55	0.46	−0.35	0.34	0.60	0.70

Table 10.

P32 equation, poor tracker with 2 tails and LT tau = 0.73, UT tau = 0.40. P25 equation, poor tracker with 2 tails and LT tau = 0.06, UT tau = 0.06.

Patient 4 poor tracker, UT R² squared (non adj, adjusted) = (0.57, 0.61).

Patient 25 – good tracker, 2 tails, R² squared (non adj, adjusted) = (0.60, 0.70).

Score = 1.55 + 0.46Dose-0.35Dose*LT + 0.34Dose*UT (Table 10): The intercept of 1.55 indicates that the patient is experiencing “chronic” background agitation. The slope of 0.46 indicates that in the moderate agitation zone nurses tend to strongly underestimate the patient’s agitation severity. In severe agitation zone, nurses tend to overestimate the agitation severity. In the mild agitation zone, nurses tend to assign a rating that is, on average, 0.35 points lower than expected for the patient’s given agitation severity. In the severe agitation zone, nurses tend to overestimate the agitation severity on average, 0.34 points higher than expected for the patient’s given agitation severity (Figure 9, see also Figure 11). For P25 LT tau = 0.06, UT tau = 0.06, and its (LT, UT) dose threshold = (2.68, 4.43), (LT, UT) score threshold of (2.41, 3.94).

Score = 0.79 + 0.52Dose+0.45Dose*UT (Table 11): The slope of 0.52 indicates that in the moderate agitation zone nurses tend to strongly underestimate the patient’sagitation severity. In the severe agitation zone, nurses tend to overestimate the agitation severity. There are 14 (out of a total of 63) such occurrences indicating that approximately only one in every four ratings, nurses tend to overestimate patients’ agitation severity (Figure 13, see also Figure 11). For P4 UT tau = 0.69, and its UT dose threshold = 2.82 and UT nurses’ score threshold = 2.75.

Patient 4 - poor, Upper Tail; Adjusted R-squared: 0.61				Estimate	Standard Error	t value	P value
Intercept				0.79	0.446	1.77	0.082		.
Dose				0.521	0.238	2.19	0.032		*
Dose*UT				0.445	0.189	2.35	0.022		*
Patient 23 - good, Upper Tail; Adjusted R-squared: 0.68			Estimate		Standard Error	t value	P value
Intercept			0.343		0.257	1.34	0.187		NS
Dose			0.576		0.24	2.4	0.0199		*
Dose*UT			0.524		0.196	2.68	0.0098		**
Patient 37 - good, Upper Tail; Adjusted R-squared: 0.75			Estimate		Standard Error	t value	P value
Intercept			0.165		0.273	0.61	0.546		NS
Dose			0.778		0.186	4.18	5.60E−05		***
Dose*UT			0.318		0.143	2.23	0.028		*
		Linear correlation (r)			Novel Regression (estimates)				adj R^2quare
P no.	Tail status	Main region	Lower tail	Upper tail	Intercept	Dose Slope	Dose Slope Lower	Dose Slope Upper	Simple LR R^2	Novel method R^2
4	Upper	0.38		0.66	0.79	0.52		0.45	0.57	0.61
23	Upper	0.51		0.62	0.34	0.58		0.52	0.64	0.68
37	Upper	0.54		0.75	0.17	0.78		0.32	0.74	0.75

Table 11.

P4 equation, poor tracker with upper tail and UT tau = 0.69. P23 equation, poor tracker with upper tail and UT tau = 0.68.

Figure 13.

P4 poor tracker with upper tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P4 nurses’ A-S score with dose. UT tau = 0.69, UT dose threshold = 2.82 and UT nurse score threshold = 2.75.

Patient 23 good tracker, UT R² squared (non adj, adjusted) = (0.64, 0.68).

Score = 0.34 + 0.58Dose+0.52Dose*UT (Table 11): The slope of 0.58 indicates that in the moderate agitation zone nurses tend to strongly underestimate the patient’s agitation severity. In periods of severe agitation, nurses tend to overestimate the agitation severity. There are 10 (out of a total of 57) such occurrences indicating that approximately only one in every five ratings, nurses tend to overestimate patient 23’s agitation (Figure 14, see also Figure 11). For P23 UT tau = 0.68, and its UT dose threshold = 1.77 and UT nurses’ score threshold = 1.97.

Figure 14.

P23 good tracker with upper tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P23 nurses’ A-S score with dose. UT tau = 0.68, UT dose threshold = 1.77 and UT nurse score threshold = 1.97.

Patient 37 good tracker, UT R² squared (non adj, adjusted) = (0.74, 0.75).

Patient 4 - poor, Upper Tail Score = 0.17+ 0.78Dose+0.32Dose*UT (Table 11): The slope of 0.78 indicates that in the moderate agitation zone nurses tend to moderately underestimate the patient’s agitation severity. In the patient’s severe agitation zone, nurses tend to overestimate agitation severity. There are 25 (out of a total of 123) such occurrences indicating that approximately only in one in every five ratings, nurses tend to overestimate patients’ agitation severity accurately. (Figure 10, see also Figure 12). For P37 UT tau = 0.77, and its UT dose threshold = 2.53 and UT nurses’ score threshold = 2.99.

P37 equation, good tracker with upper tail and UT tau = 0.77.

Patient 2 poor tracker, LT R² squared (non adj, adjusted) = (0.42, 0.44).

Score = 0.27 + 0.85Dose - 0.37Dose*LT (Table 12): The slope of 0.85 indicates that in the moderate agitation zone nurses tend to underestimate the patient’s agitation severity. In the mild agitation zone, nurses tend to assign a rating that is, on average, 0.37 points lower than expected for the patient’s given agitation severity. There are 28 out of a total of 63 such occurrences indicating that in around one in every two ratings, nurses tend to underestimate agitation severity (Figure 15, also Figure 11).

Patient 2 - poor, Lower Tail; Adjusted R-squared: 0.44			Estimate		Standard Error	t value		P value
Intercept			0.269		0.208	1.29		0.201	NS
Dose			0.85		0.158	5.36		1.40E−06	***
Dose*LT			−0.367		0.184	−1.99		0.052	.
Patient 15 - good, Lower Tail; Adjusted R-squared: 0.47			Estimate		Standard Error	t value		P value
Intercept			0.861		0.962	0.9		0.37412	NS
Dose			0.834		0.206	4.05		0.00015	***
Dose*LT			−0.344		0.177	−1.94		0.05664	.
		Linear correlation (r)			Novel Regression (estimates)				adj R^2quare
P no.	Tail status	Main region	Lower tail	Upper tail	Intercept	Dose Slope	Dose Slope Lower	Dose Slope Upper	Simple LR R^2	Novel method R^2
2	Lower	0.39	0.77		0.27	0.85	−0.37		0.42	0.44
15	Lower	0.28	0.75		0.86	0.83	−0.34		0.45	0.47

Table 12.

P2 equation, poor tracker with lower tail and LT tau = 0.64. P15 equation, good tracker with lower tail and LT tau = 0.74.

Figure 15.

P2 poor tracker with lower tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P2 nurses’ A-S score with dose. LT tau = 0.64, with LT dose threshold = 0.96 and LT nurses’ score threshold = 0.83.

Patient 15 good tracker, LT R² squared (non adj, adjusted) = (0.45, 0.47).

Score = 0.86 + 0.83Dose - 0.34Dose*LT (Table 12): The slope of 0.83 indicates that in the moderate agitation zone nurses tend to underestimate the patient’s agitation severity accurately. In the mild agitation zone, nurses tend to underestimate the patient’s agitation severity even more. There are 22 (out of a total of 43) such occurrences indicating that approximately one in every three ratings, nurses tend to strongly underestimate the patients’ agitation severity (Figure 16, see also Figure 11).

Figure 16.

P15 good tracker with lower tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P15 nurses’ A-S score with dose. LT tail tau = 0.74, LT dose threshold = 3.47 and LT score threshold = 3.05.

Patient 28 poor tracker, LT R² squared (non adj, adjusted) = (0.57, 0.59) Score = 0.01 + 0.97Dose-0.35Dose*LT (Table 13).

Patient 28 – poor, Lower Tail			Estimate		Standard Error	t-value	P-value		adj R^2 = 0.64
Intercept			0.00876		0.17549	0.05	0.96023
Dose			0.96582		0.07356	13.13	< 2e−16		***
Dose*LT			−0.345		0.1014	−3.4	0.00081		***
Patient 35 – poor, Lower Tail			Estimate Std.		Error t	t value	P value		adj R^2 = 0.644
Intercept			0.2322		0.1128	2.06	0.04084		*
Dose			0.9097		0.0597	15.24	< 2e−16		***
Dose*LT			−0.7898		0.2253	−3.51	0.00056		***
Patient 27 – poor, no tails			Estimate Std.		Error t	t value	P value		adj R^2 = 0.71
Intercept			−0.3644		0.2062	−1.77	0.079		.
Dose			1.061		0.0453	23.4	<2e-16		***
		Linear correlation (r)			Novel Regression (estimates)				adj R^2 quare
P no.	Tail status	Main region	Lower tail	Upper tail	Intercept	Dose Slope	Dose Slope Lower	Dose Slope Upper	Simple LR R^2	Novel method R^2
28	Lower	0.46	0.84		0.00	0.97	−0.35		0.57	0.59
35	Lower	0.72	0.77		0.23	0.91	−0.79		0.63	0.64
27	None	0.84			−0.36	1.06			0.71

Table 13.

P28 equation, poor tracker with lower tail LT tau = 0.79. P35 equation, poor tracker with lower tail LT tau = 0.74. P27 equation, poor tracker no tails., tau = 0.78 (Figures 18 and 19).

The slope of 0.97 indicates that in the moderate agitation zone nurses tend to estimate the patient’s agitation severity fairly accurately. In the mild agitation zone, however, nurses tend to underestimate patient 28’s agitation severity more. There are 29 (out of a total of 203) such occurrences indicating that one in every seven ratings, nurses tend to strongly underestimate patients’ agitation severity (Figure 17, see also Figure 18).

Figure 17.

Line plots of nurses’ score (observed, red) vs. dose (black line), with 95% WPB bands for P27 and P35 poor trackers with and P27 poor tracker with no tails.

Figure 18.

Line plots of nurses’ score (observed, red) vs. dose (black line), with 95% WPB bands for P28 and P35 both poor trackers with lower tails.

Figure 19.

P35 poor tracker with a lower tail - copula plot, main tau, K- plot and tail tau (upper panel); best fit line(s) (lower panel) relating P35 nurses’ score with dose. LT tail tau = 0.74, LT dose threshold = 0.72 and LT score threshold = 0.66.

Patient 35 poor tracker, LT R² squared (non adj, adjusted) = (0.63, 0.64) Score = 0.23 + 0.91Dose-0.79Dose*LT (Table 13): The intercept of 0.23 indicates that the patient is experiencing mild “chronic” background agitation. The slope of 0.91 indicates that in the moderate agitation zone nurses tend to slightly underestimate the patient’s agitation severity. In the mild agitation zone, nurses tend to strongly underestimate the patient’s agitation severity even more. There are 42 (out of a total of 211) such occurrences indicating that approximately one in every five ratings, nurses tend to underestimate patients’ agitation severity (Figures 8 and 18).

4. Conclusions

Copulas were successfully implemented to capture non-linear dependence, and establish the presence of lower (LT) and/or upper tail (UT) dependence between the nurses’ A-S rating and the automated sedation dose (Figure 20). Establishing the presence of tail dependence and patient-specific lower and/or upper thresholds for areas with different agitation intensities has significant implications for effective administration of sedatives. Copulas unique to the poor trackers were the Survival BB7 and Rotated Tawn type 1, 180. Unique copula types of the good trackers were Clayton, BB8, and t copula. To the best fit copulas, we established for each patients’ bivariate score and dose trajectories regions of mild, moderate, and severe agitation and their lower and upper tail thresholds, if any, in the dependence, relationship via K-plots and our novel equation relating patient’s nurses’ score to dose. We found that for lower tail dependence, nurses tend to underestimate the patient’s agitation in the moderate agitation zone. In the mild agitation zone, nurses tend to assign a rating that is, on average, 0.30 to 0.45 points lower than expected for the patient’s given agitation. For upper tail dependence nurses tend to either moderately or strongly underestimate patient’s agitation in the moderate agitation zone; but in periods of severe agitation, tend to overestimate a patient’s agitation. When both lower and upper tails exist, nurses tend to strongly underestimate agitation in the moderate zone and then tend to still underestimate agitation in mild agitation periods; but in the severe zone nurses tend to overestimate agitation on average by 0.34 points. We also determined the number of occurrences over the patient’s total ICU stay, where nurses tend to under/over- estimate agitation. Finding thresholds and regions of mismatch between the nurses’ scores and sedation dose potentially provides a way to improved alerting systems for over/under-sedation. Our approach lends credence to augmenting conventional RASS and SAS agitation measures with semi-automated systems and identifying thresholds and regions of deviance for alerting increased risk.

Figure 20.

P27 poor tracker with no tails- copula plot, main tau = 0.70, K- plot and tail tau (upper panel); best fit line (lower panel) relating P27 nurses’ score with dose.

			Lower Tail tau (τ)		Upper Tail tau (τ)	Dose threshold			Nurses’ Score threshold
P no.	Tails	n	n	Per centile	n	Per centile	Lower	Upper	Lower	Upper
1	Lower	59	19	32.20			0.52		0.38
2	Lower	63	28	44.44			0.96		0.83
3	Lower	63		22.22			1.35		1.31
4	Upper	63			15	76.19		2.82		2.75
5	Lower	48	17	35.40			1.25		1.06
6	Lower	63	23	36.51			1.38		1.23
7	Lower	63	8	12.70			0.44		0.15
8	Lower	127	33	25.98			3.02		2.57
10		255
11		111
12		127
13		63
14	Lower	49	20	40.82			0.96		0.82
15	Lower	63	22	34.92			3.47		3.05
16		127
17		63
18		63
19	Lower	127	28	22.04			6.46		5.23
20		127
21		61
22		127
23	Upper	57			10	82.46		1.77		1.97
24		127
25	Both	63	14	33.33	15	68.25	2.68	4.43	2.41	3.94
26		63
27		223
28	Lower	203	59	29.06			1.70		1.20
29		53
30		60
31		255
32	Both	252	44	20.64	38	74.60	3.90	7.02	3.62	7.63
33		255
34		127
35	Lower	211	42	19.91			0.72		0.66
36		57
37	Upper	123			17	86.18		2.53		2.99

Table A.1.

List of copulas and upper and/or lower tail positions (percentiles) and associated lower/upper thresholds for dose and nurses’ A-S score.

Shaded rows are the poor trackers.

		Linear correlation (r)			Novel Regression (estimates)				adj R^2
P no.	Tail status	Main region	Lower tail	Upper tail	Intercept	Dose Slope	Dose Slope Lower	Dose Slope Upper	Simple LR R^2	Novel method R^2
1	Lower	0.44	0.81		0.152	0.84	−0.46		0.49	0.53
2	Lower	0.39	0.77		0.27	0.85	−0.37		0.42	0.44
3	Lower	0.36	0.65		0.91	0.58	−0.50		0.30	0.34
4	Upper	0.38		0.66	0.79	0.52		0.45	0.57	0.61
5	Lower	0.21	0.75		0.49	0.75	−0.45		0.33	0.39
6	Lower	0.37	0.56		0.66	0.68	−0.36		0.36	0.41
7	Lower	0.85	0.99		−0.22	1.17	−0.45		0.77	0.77
8	Lower	0.64	0.76		0.01	1.01	−0.29		0.62	0.63
10		0.82			−0.51	1.19			0.67
11		0.83			−0.19	1.21			0.68
12		0.68			0.18	0.94			0.46
13		0.87			−0.43	1.12			0.75
14	Lower	0.50	0.70		0.11	0.95	−0.40		0.61	0.64
15	Lower	0.28	0.75		0.86	0.83	−0.34		0.45	0.47
16		0.81			−0.17	1.03			0.66
17		0.71			0.19	0.90			0.50
18		0.77			0.01	1.06			0.58
19	Lower	0.47	0.55		0.80	0.94	−0.34		0.48	0.50
20		0.70			−0.31	1.16			0.49
21		0.74			−0.10	1.08			0.55
22		0.77			−0.58	1.16			0.59
23	Upper	0.51		0.62	0.34	0.58		0.52	0.64	0.68
24		0.78			−0.50	1.14			0.60
25	Both	0.44	0.41	0.00	1.55	0.46	−0.35	0.34	0.60	0.70
26		0.78			−0.07	1.08			0.61
27		0.84			−0.36	1.06			0.71
28	Lower	0.46	0.84		0.008	0.97	−0.35		0.57	0.59
29		0.86			0.24	0.89			0.73
30		0.55			−0.04	1.04			0.29
31		0.79			1.00	1.16			0.63
32	Both	0.26	0.55	0.61	3.23	0.43	−0.75	0.33	0.60	0.67
33		0.83			−0.98	1.23			0.68
34		0.74			−9.23	1.06			0.55
35	Lower	0.72	0.77		0.23	0.91	−0.79		0.63	0.64
36		0.86			−0.45	1.11		0.28	0.73
37	Upper	0.54		0.75	0.17	0.78		0.32	0.74	0.75

Table A.2.

Upper and/or lower tail dose relationships, equations, and associated change in adjusted R² for simple LR vs our novel approach.

Shaded rows indicate the poor trackers. Boxes regarding indicate largest differences between R² equation fit.