|Home | About | Journals | Submit | Contact Us | Français|
Controversy has attended the relationship between risk-adjusted mortality and process-of-care. There would be advantage in the establishment, at the data-base level, of global quantitative indices subsuming the diversity of process-of-care.
A retrospective, cohort study of patients identified in the Australian and New Zealand Intensive Care Society Adult Patient Database, 1993-2003, at the level of geographic and ICU-level descriptors (n = 35), for both hospital survivors and non-survivors. Process-of-care indices were established by analysis of: (i) the smoothed time-hazard curve of individual patient discharge and determined by pharmaco-kinetic methods as area under the hazard-curve (AUC), reflecting the integrated experience of the discharge process, and time-to-peak-hazard (TMAX, in days), reflecting the time to maximum rate of hospital discharge; and (ii) individual patient ability to optimize output (as length-of-stay) for recorded data-base physiological inputs; estimated as a technical production-efficiency (TE, scaled [0,(maximum)1]), via the econometric technique of stochastic frontier analysis. For each descriptor, multivariate correlation-relationships between indices and summed mortality probability were determined.
The data-set consisted of 223129 patients from 99 ICUs with mean (SD) age and APACHE III score of 59.2(18.9) years and 52.7(30.6) respectively; 41.7% were female and 45.7% were mechanically ventilated within the first 24 hours post-admission. For survivors, AUC was maximal in rural and for-profit ICUs, whereas TMAX (≥ 7.8 days) and TE (≥ 0.74) were maximal in tertiary-ICUs. For non-survivors, AUC was maximal in tertiary-ICUs, but TMAX (≥ 4.2 days) and TE (≥ 0.69) were maximal in for-profit ICUs. Across descriptors, significant differences in indices were demonstrated (analysis-of-variance, P ≤ 0.0001). Total explained variance, for survivors (0.89) and non-survivors (0.89), was maximized by combinations of indices demonstrating a low correlation with mortality probability.
Global indices reflecting process of care may be formally established at the level of national patient data-bases. These indices appear orthogonal to mortality outcome.
The outcomes paradigm is now a dominant influence within medicine  and critical care is no exception to this movement . The 1986 paper by Knaus et al , evaluating outcomes of a cohort of 13 intensive care units (ICU), established the notion of institutional or provider performance within the critical care discipline by way of the nexus between risk-adjusted mortality and process-of-care, the latter established through questionnaire, on-site visit and case-note review. Similar investigations were concurrently reported in the general medical literature by Dubois and colleagues . A discordant debate has subsequently occurred regarding the relationship between risk-adjusted mortality and process-of-care, the latter being variously assessed [5,6]. On the one hand mortality "...is unlikely to be a sufficient statistic for quality" ; yet, the felicity with which process may be measured is no guarantee that "measuring ...process and reporting performance will improve outcomes" . Contra-wise to the relationship of process-of-care and mortality outcome, a recent study has suggested that the "notion that hospitals with higher risk-adjusted mortality rates have poorer quality care is neither consistent nor reliable" . However, there is a certain circularity in these arguments: reliance on outcome measures (mortality or length of stay ) is criticised from the stand point of process-of-care  (adherence to checklists ), which in turn finds its (ultimate) assessment in terms of the effect on precisely those outcomes which have been rejected in the first place.
As opposed to a piecemeal examination of single indicators or a composite-scores approach , there would appear to be advantage in the establishment of global quantitative indices [14,15] which would subsume the diversity of process-of-care and avoid the necessity of direct examination of the modalities of the latter . We sought to establish such indices at the level of a bi-national intensive care patient data-base [16,17], the Adult Patient Database (APD) of the Australian and New Zealand Intensive Care Society (ANZICS); by analysis of:
• the hazard of patient hospital discharge, estimated using time-to-event analysis , as reflecting the time-course of process-of-care. The components of the time-hazard curve were determined using pharmaco-kinetic methods.
• individual patient ability to maximize output , in this case length of stay, for a given set of physiological inputs, the individual patient component variables of the Acute Physiology and Chronic Health Evaluation (APACHE) III severity of illness algorithm . This ability was conceptualised as one of technical production efficiency ("economic" efficiency output/input, scaled [0,1] ), estimated by the econometric technique of stochastic frontier analysis (SFA) .
We also determined the degree of correlation, or independence (orthogonality), between these global process-of-care indices and mortality.
As previously described , the ANZICS APD  was interrogated to define an appropriate patient set over the time period 1993-2003. In brief, physiological variables collected, in accordance with the requirements of the APACHE III algorithm , were the worst in the first 24 hours after ICU admission, and all first ICU admissions to a particular hospital for the period 1993-2003 were selected. Records were used only when all three components of the Glasgow Coma Score (GCS) were provided; records for which all physiologic variables were missing were excluded, and for the remaining records, missing variables were replaced with the normal range and weighted accordingly. ICU and hospital length-of-stay, initially recorded in hours, were transformed to fractional days. Patients with an ICU length of stay > 60 days and hospital length of stay > 365 days were not considered in formal analysis. Exclusions: unknown hospital vital outcome and date of discharge; patients with an ICU length of stay ≤ 4 hours; and patients aged < 16 years of age. Access to the data was granted by the ANZICS Database Management Committee in accordance with standing protocols; local hospital (The Queen Elizabeth Hospital) Ethics of Research Committee approval was waived.
ICU-year units were formed by ICU-site × calendar-year interaction with a minimum patient number set at 150 to ensure estimation stability. Categorical predictors were parameterized as indicator variables with the reference level (= 0) indicated in parentheses in the following list:
• mechanical ventilation (not ventilated)
• ICU-level, as defined in the ANZICS database, as Rural, Metropolitan, Tertiary and Private (Tertiary)
• Geographical-location; that is New Zealand and the States of the Commonwealth of Australia, excluding Western Australia (New South Wales (NSW), the largest contributor). A reference population-density map is provided in Additional file 1.
• patient surgical status as post-elective surgery, post-emergency surgery and non-surgical (non-surgical)
• annual patient admission number (reflecting ICU "size"), created by the ICU-site × calendar-year interaction, was categorized at the median, the reference category being that denoting the highest number of yearly admissions
Variables were reported as mean (SD), except where otherwise indicated. Distributional form of variables of interest was displayed using kernel density estimates and the empirical distribution and parameter 95% confidence intervals (CI) were computed via bootstrapping (1000 repetitions) . In the presence of skewed distribution, median point-estimates were reported. Group differences between continuous variables were estimated using parametric or non-parametric analysis of variance (ANOVA, adjusted (Sidak) for multiple comparisons) where indicated. Stata™ (Version 10 MP, 2007; College Station, Texas) statistical software was used.
Hospital mortality probabilities (MP), summed over descriptors, were estimated using a two-level, patients within ICU-year units, random effects [random-intercept] logistic regression model as previously reported . Descriptor standardised mortality rates (SMR) and 95%CI were estimated using the parametric approach of Rapoport et al .
The unit of analysis was patients within descriptors. For hospital survivors, length of stay of non-survivors was defined as >> the maximum length-of-stay of alive discharges, with no censoring. The hazard of hospital (or ICU) discharge was computed via time-to-event analysis with kernel density smoothing . Similarly, for hospital non-survivors, length of stay of survivors was defined >> maximum length-of-stay of those dying, with no censoring . Such an approach obviated the analytical problem of competing risks . The smoothed hazard estimates and 95%CI, truncated at 30 days and adjusted for the mean values of the covariates age and APACHE III score for each descriptor, were output as individual data files and separately processed using standard pharmaco-kinetic techniques (the Stata™ "pk" suite of commands ) to estimate parameters of the "hazard" profile; in particular: peak hazard (CMAX) and time to peak hazard (TMAX), reflecting the initial intensity of and time to the maximum (rate of) patient hospital discharge, respectively; area under the hazard curve (AUC), reflecting the integrated experience of the discharge process; and the "elimination rate" (KE), the negative of the parameter estimate for a linear regression of log(time) on hazard, which determines the half-life (t1/2) of the overall hazard-profile (). The justification for this approach was that an initial (random effects) first-order compartment model  provided a good fit to the data (see Additional file 2). These parameter estimates were understood as global indices reflecting aspects of descriptor process-of-care.
Patients were considered as producers, seeking to avoid waste by obtaining maximised outputs from given inputs or, by minimizing input use in the production of given outputs ; where, in this context, "maximize" is used in the sense of optimize. The notion of productive efficiency corresponds to technical efficiency, the latter being estimated by SFA. A stochastic production frontier model may be estimated, using the Stata™ module "frontier", as a log-linear function (f): ; where TEi = exp(-ui) and ui > 0, here assumed exponentially distributed; yi is ICU/hospital length of stay; xijs are (logged) acute physiology score and chronic health evaluation variables (age, GCS, temperature, heart rate, arterial pH, arterial PaCO2, creatinine, mean arterial pressure, white cell count, plasma bilirubin, plasma glucose and total (APACHE III) co-morbidity count); , i = 1, ..., N (idiosyncratic patient component). Adjustment for heteroscedasticity of the variance function of both u and v (as provided for by the Stata™ module "frontier") was undertaken in model development with a combination of appropriate patient (gender, patient surgical status), treatment (ventilation status) and provider descriptor variables (calendar year, ICU level, annual patient admission number and geographical locality). Patient efficiencies (TE, scaled [0, 1], where 1 the optimal production frontier) were summed over categorical descriptors of interest.
The multivariate relationships (joint distribution) between the indices (MP, TE, TMAX, AUC, CMAX and KE), for survivors and non-survivors across descriptors, were displayed using biplots . Biplots consist of lines, reflecting the dataset variables, and "dots" to show the observations. The length of the lines approximates the variances of the variables (the longer the line, the higher is the variance) and the cosine of the angle between the lines approximates the correlation; the closer the angle is to 90, or 270 degrees, the smaller the correlation (orthogonality); an angle of 0 or 180 degrees reflecting a correlation of 1 or -1, respectively. Variable inclusion in the biplots was adjusted to maximize the explained variance.
The data set consisted of 223129 patients from 99 ICUs over an 11 year period. Mean (SD) age and APACHE III score were 59.2 (18.9) years and 52.7 (30.6) respectively; 41.7% were female and 45.7% were mechanically ventilated within the first 24 hours post ICU-admission. Overall ICU and hospital mortalities were 10.4% and 16.1% respectively. ICU length of stay was 3.6 (5.6) (median 1.8, interquartile range 2.9 [0.9-3.8]) days and hospital length of stay was 16.4 (19.5) (median 10.1, interquartile range 14.6 [5.1-19.7]) days. Patient categorization was non-operative (55.2%), elective surgical (28.7%) and emergency surgical (16.1%). Annual patient admission number created by the ICU-site × calendar-year interaction, was categorized at the median (711 patients); this categorization was further used to create "ICU-hospital/geographical area/size" descriptors (n = 35).
Hazard-of-discharge curves, for geographical areas and hospital-ICU levels, are seen for alive-discharges and those dying in Figures Figures11 and and22 respectively. The curves for those dying were not unexpectedly displaced to the left compared with those surviving. Initial analysis of the smoothed hazard curves used a non-linear mixed effects approach which yielded parameter estimates displaying good between-descriptor discriminative properties (Additional file 2). As this analysis required specialised software , standard pharmacokinetic parameters, with bootstrapped 95% CI, were calculated. Descriptor technical efficiencies (Figure (Figure3)3) and mortality probabilities revealed non-normality and median point estimates were reported. Complete summaries at the level of geographical × hospital-ICU level categories × ICU "size", with standardised mortality rates, are given for both survivors and non-survivors in Tables Tables11 and and22 respectively.
Summary statistics for process-of-care indices and mortality probabilities, for both survivors and non-survivors, are displayed for descriptors ICU "size" (> 711, < 711 admissions per year); ICU level by "size" and geographical-location by "size" in Tables Tables33 and and4.4. For both survivors and non-survivors, indices demonstrated statistically significant differences between the individual categories of various descriptors: ICU "size" (> 711, < 711 admissions per year); ICU level; ICU level by "size"; geographical-location by "size"; and geographical × hospital-ICU level categories × ICU "size" (ANOVA, P ≤ 0.0001); albeit these differences reflected the large size of both the initial data-set and the empirical distributions of the indices (1000 bootstrap samples). Maximal values of indices were seen (Tables (Tables33 and and4;4; ANOVA, P < 0.001 compared with all other ICU levels) for:
the descriptor ICU level by "size":
1. for survivors, at the rural and private level (both < 711 and >711 yearly admissions) for AUC; at tertiary levels (both < 711 and >711 yearly admissions) for TMAX and TE; at rural (>711 yearly admissions) and metropolitan (<711 yearly admissions) levels for KE; at rural and private levels for CMAX.
2. for non-survivors, at the tertiary level (both < 711 and >711 yearly admissions) for AUC; at the private level for TMAX; at rural (<711 yearly admissions) and metropolitan (<711 yearly admissions) levels for KE; at the tertiary levels (both < 711 and >711 yearly admissions) for CMAX
the descriptor ICU "size":
1. for survivors: for ICUs with > 711 admissions per year for AUC and TE; for ICUs with < 711 admissions per year for TMAX
2. for non-survivors: for ICUs with > 711 admissions per year for AUC; for ICUs with < 711 admissions per year for TMAX and TE.
The relationships between the global indices are demonstrated in Figure Figure4.4. For survivors (left panel), the total explained variance was 0.886, and the relationships appeared discrete: both technical efficiency (TE) and time to maximal hazard (TMAX) were coincident and tended to be orthogonal to indices reflecting the maximum hazard of alive (hospital) discharge (CMAX) and the total discharge experience of the patient (AUC). When included in the plot, mortality probability (and SMR) was, not surprisingly, directionally opposite to that of AUC and did not increase the explained variance. For non-survivors (right panel), the total explained variance was 0.892. AUC and mortality probability were almost coincident and directionally equivalent to CMAX, but orthogonal to both TMAX and TE. The SMR and mortality probability were coincident (not shown). The "elimination rate" (KE) added no increment to the total explained variation for survivors or non-survivors.
We have proffered a number of global indices, their uncertainty estimates and multivariate relationships, as reflecting the integrated discharge process for both survivors and non-survivors at the level of various data-base descriptors. For each of the indices, we were able to establish formal statistical difference between, and characteristic clustering within, database descriptor categories for both survivors and non-survivors. These indices would appear to be both internally consistent and plausible proxies for the underlying process(es)-of-care which, in concert with patient characteristics, determine the shape of the particular descriptor hazard of discharge.
Our motivating concepts were two-fold:
1. that of "conditional length of stay", introduced by Silber and colleagues, as the length of stay after a stay is prolonged, reflecting both patient complications and/or co-morbid medical conditions, and provider ability to manage complicated cases [30,31]. As suggested by Silber et al, and based upon empirical analyses: "By studying the CLOS, one can determine when the rate of hospital discharge begins to diminish-without the need to directly observe complications. Policy makers looking for an objective outcome measure may find that CLOS aids in the analysis of a hospital's management of complicated patients..." . The analytic focus was restricted to the prolongation point or day, as estimated by the Hollander and Proschan statistic . Our development of this approach involved modelling the hazard-profile as a "concentration" curve with the estimation of standard pharmaco-kinetic based parameters which would allow a more complete description of the hazard-time curve, in the sense of modelling underlying processes, and allowing formal statistical comparisons.
2. the notion of efficiency. In an innovative study of patients with severe head trauma, Nathanson et al  used data envelopment analysis (DEA) to calculate individual patient "efficiency " scores ([0 to 1]) based upon the ability to maximize output (in this case, cerebral perfusion pressure) for a given set of physiological inputs. Patients with high efficiency scores were found to have an improved functional outcome on ICU discharge. At the individual ICU level, Junoy  used DEA to establish an efficiency frontier across the bivariate relationship between severity adjusted survival (effectiveness or output) and length of stay (resources or input) and to compute technical efficiency of "output quality". Length of stay, raw, adjusted or observed minus expected, has been used as an indicator of ICU/hospital performance, from both an economic  and quality  perspective. Traditionally, length of stay has been estimated by ordinary least squares regression (or its variants) by maximising the mean expectation; the residual (difference between observed and predicted) being interpreted, in the current context, as arising from inefficiency. We applied the concept of "efficiency" to individual patient length of stay to model technical production efficiency using a parametric SFA approach. The latter separates the residual into an inefficiency component (ui, positive departures from the (best practice) production frontier) and all other sources of model error (vi), such as random shocks and measurement error [15,36]. To this extent it is less sensitive to outliers than DEA, a deterministic non-parametric technique, assuming no measurement error and requiring a more rigid best-practice production frontier based upon a small subset of efficient peers .
That these indices captured aspects of descriptor process-of-care was suggested by the analysis of maximum values. For instance, tertiary level ICUs where more complex and severely ill patients were located  and human and material resources were presumably greatest, demonstrated maximum TE for survivors and AUC for non-survivors; whereas, at the rural and for-profit ICU levels, where such patient and resource conditions did not necessarily obtain, AUC was maximal for survivors. With respect to ICU "size", maximal TE was located in tertiary ICUs (> 711 yearly admissions) for survivors; but of interest, for non-survivors, was located in Private ICUs (< 711 admissions per year). These "size" effects must be interpreted against the favourable effect of ICU "size" < 711 yearly admissions on hospital mortality (OR 0.84, P < 0.0001 compared with > 711 yearly admissions) adduced in a previous analysis of the same data-base . The regional-geographic differences in these indices (Tables (Tables1,1, ,2,2, ,3,3, ,4,4, above) presumably reflect both the particular nature of our data-base and determinants such as the distribution of human [37,38] and non-human resources  and socio-economic factors.
Analysis of the joint distribution of the global process-of-care indices revealed a weak correlation with, or orthogonality to, mortality outcome (SMR and/or mortality probability). This would appear to be the first formal demonstration of such a (lack-of) relationship and is consistent with the literature reviews of Thomas and Hofer  and Pitches et al , although both reviews noted considerable and possibly confounding between-study heterogeneity in recorded process-of-care and mortality risk-adjustment measures. The low sensitivity of individual process-of-care measures has also been previously noted . This being said, any such formal "independence" of these measures would have important implications for both health policy and the design and interpretation of trials assessing process-of-care interventions .
Our postulates are predicated upon the utility of both the pharmaco-kinetic and efficiency analyses. The former would appear to have basis in the good fit of a first-order compartment model; this approach may be formally extended to embrace non-linear mixed effects modelling . Controversy has attended the appropriate form of efficiency analysis  and the use of either DEA or SFA; much of this criticism is context dependent; for instance, the efficiency of public services where multiple outputs occur. Our use of SFA was motivated by modelling flexibility and potential extensions to accommodate random effects and different distributions for the inefficiency component [44,45].
Our analysis, at the level of a bi-national database, lacks the empirical grounding of the investigations of Silber and co-workers [30,31] mentioned above and would require such validation. However, at the analytic level, these indices appeared both consistent and intuitively reasonable. In the interests of parsimony the primary descriptors were geographical and hospital-ICU level categories and the data-set was considered as a single cross-sectional unit; that is no analysis by calendar year was undertaken, although this ensured estimation stability. Similarly, we considered only hospital outcomes and truncated analysis time to 30 days, the latter again for estimation stability in the output files from the smoothed hazard estimates. However, our approach is rich with possibilities for various extensions to the individual ICU level and ICU hazard-of-discharge, longitudinal analysis , and obviously, different non-ICU patient cohorts.
Global indices reflecting process of care may be formally established at the level of national patient data-bases, thus allowing comparisons between providers/descriptors. These indices appear orthogonal to mortality outcome; such a relationship would have implications for health care policy and the design and interpretation of trials assessing process-of-care interventions.
The authors declare that they have no competing interests.
The study was conceived, designed, (data)-analysed, written and critically revised jointly by both authors (JLM, PJS).
The pre-publication history for this paper can be accessed here:
Reference population-density map of Australia and New Zealand. Population density map (inhabitants per square kilometre) of Australia (left figure) and New Zealand (right figure); not drawn to scale. Approximate land masses: Australia, 7706143 km2; New Zealand 270626 km2. Gray-scale legends indicate inhabitant number (range) per km2. WA, Western Australia; NY, Northern Territory; QLD, Queensland; NSW, New South Wales; VIC, Victoria; TAS, Tasmania.
Random effects first-order compartment model. A brief explanation of the self-starting first-order compartment model. Figure: Parameter estimates with 95%CI for a first-order compartment non-linear (self-starting "ssfol") mixed effects model for geographical descriptors. TAS, Tasmania; NSW, New South Wales; QLD, Queensland; VIC, Victoria; NT, Northern territory; NZ, New Zealand.