Risk factors and action thresholds for the novel coronavirus pandemic. Insights from the Italian Society of Nephrology COVID-19 Survey.

BACKGROUND AND AIM: Over 80% (365/454) of the nation's centers participated in the Italian Society of Nephrology COVID-19 Survey. Out of 60,441 surveyed patients, 1368 were infected as of April 23rd, 2020. However, center-specific proportions showed substantial heterogeneity. We therefore undertook new analyses to identify explanatory factors, contextual effects, and decision rules for infection containment.
METHODS: We investigated fixed factors and contextual effects by multilevel modeling. Classification and Regression Tree (CART) analysis was used to develop decision rules.
RESULTS: Increased positivity among hemodialysis patients was predicted by center location [incidence rate ratio (IRR) 1.34, 95% confidence interval (CI) 1.20-1.51], positive healthcare workers (IRR 1.09, 95% CI 1.02-1.17), test-all policy (IRR 5.94, 95% CI 3.36-10.45), and infected proportion in the general population (IRR 1.002, 95% CI 1.001-1.003) (all p < 0.01). Conversely, lockdown duration exerted a protective effect (IRR 0.95, 95% CI 0.94-0.98) (p < 0.01). The province-contextual effects accounted for 10% of the total variability. Predictive factors for peritoneal dialysis and transplant cases were center location and infected proportion in the general population. Using recursive partitioning, we identified decision thresholds at general population incidence ≥ 229 per 100,000 and at ≥ 3 positive healthcare workers.
CONCLUSIONS: Beyond fixed risk factors, shared with the general population, the increased and heterogeneous proportion of positive patients is related to the center's testing policy, the number of positive patients and healthcare workers, and to contextual effects at the province level. Nephrology centers may adopt simple decision rules to strengthen containment measures timely.

Introduction

The spread and dynamics of the novel coronavirus disease (COVID-19) in Italy posed serious threats to the ability and preparedness of nephrology centers to respond to the needs of patients on renal replacement therapy (RRT) during the exponential phase of the pandemic. At that time, sparse and scattered data were available on RRT patients and centers. Among the actions taken to provide clinicians and health authorities with nationwide data, the Italian Society of Nephrology launched a survey to evaluate the spread of SARS-CoV-2 in the Nephrology and Dialysis Units [1]. The study aimed to estimate the cumulative incidence of SARS-CoV-2 positive cases among hemodialysis (HD), peritoneal dialysis (PD), and renal transplant (TX) patients. It revealed some intriguing findings: higher rates of positive cases in HD than in PD and TX patients, and heterogeneous spread among Italian regions. The proportion of deceased SARS-CoV-2 positive individuals was shockingly high in HD patients, and such an unfavorable outcome also occurred more frequently in PD and TX patients than in the general population [2, 3].

The pandemic struck RRT centers and patients dramatically, and although specific protocols were soon developed and endorsed [4, 5], the absence of effective and evidence-based treatments left Nephrologists substantially defenseless. Thus, prevention and meticulous infection control remained, and perhaps still remains, the only broadly applicable barrier against SARS-CoV-2 spread and its deadly consequences.

The findings of our first article [1] were descriptive in nature and elicited three further and newer research issues that we address in the present report. First, which are the possible explanatory factors for the heterogeneous proportions of SARS-CoV-2 positivity in Italian nephrology centers? Second, to what extent is the higher frequency of SARS-CoV-2 positivity in HD patients associated with the spread of the virus within the dialysis centers? Third, how can we translate the survey findings into easily and widely applicable clinical decision rules?

To address these specific questions, we designed a “de-novo” analysis cycle based on multilevel modeling and classification trees. We used mixed-effect models to identify the explanatory factors of the heterogeneous infection spread while accounting for the control measures and policies currently being adopted by each center. Using classification and regression trees we sought to develop a simple decision algorithm to stratify the probability of SARS-CoV-2 positivity in each center.

Methods

Study oversight

The Italian Society of Nephrology (SIN) COVID-19 research group designed a nationwide survey to evaluate the impact of the COVID-19 pandemic on nephrology centers and patients. Details of the survey design, instruments, and procedures can be found in Reference 1. Briefly, on average, within 11 days, the vast majority [365/454 (80.4%)] of invited centers completed and returned the survey questionnaire, an instrument designed to obtain 17 key pieces of information about patients, workforce, and facilities during the exponential phase of the COVID-19 pandemic. All returned files were checked for consistency and merged in the survey database.

Statistical analysis

A preliminary exploratory analysis was performed to detect zero inflation in the count of SARS-CoV-2 positive patients and subsequent deaths. These latter outcomes were the model dependent variables, and the center prevalent patients represented the exposed population. We then assessed the multivariable association with the count SARS-CoV-2 positive cases among healthcare workers, and the HD cases admitted to the intensive care unit (ICU) as continuous explanatory variables. Center testing policies for patients and staff, telephone triage, in-person triage, and use of personal protective equipment (PPE) entered the models as binary co-variables. Additionally, for each center, we derived the number of days elapsed since the local lockdown date [6], and the cumulative rate of positive cases in the local general population as of April 23rd [7]. Since centers are clustered in provinces and provinces in regions, three multilevel models were considered. Model 1 included provinces as the second level, model 2 included regions as the second level, and model 3 provinces nested into regions. We used both the multilevel Poisson regression and the multilevel negative binomial regression for each model, and we chose the best model using the likelihood-ratio (LR) test, when applicable, and the Akaike’s information criterion (AIC).

For the multilevel model, we estimated the variance partition coefficient (VPC) that is the proportion of variation that is beyond what explained by the fixed predictors and informs on the existence of a contextual effect [8]. Fixed predictor coefficients were exponentiated and reported as cumulative incidence rate ratios (IRR) with 95% confidence intervals (CI).

We used Classification and Regression Tree (CART) analysis [9-11] to develop a simple and interpretable set of rules to support clinical decision making. Tree-based models are a class of nonparametric algorithms that work by partitioning the feature space into a number of non-overlapping regions with similar response values using a set of splitting rules.

The presence or absence of SARS-CoV-2 positivity was the target variable. A tenfold cross-validation method was used to evaluate the model reliability. Statistical analyses were performed using STATA 15 (StataCorp LLC, TX, USA), and R Software Ver 4.0.2 (R Core Team, 2020).

Geospatial analysis

Geographical Information System (GIS) technology was used to produce multilevel maps of the model predictions on Italian provinces using the Istituto Nazionale di Statistica (ISTAT) shapefiles [12]. Geospatial mapping was performed using R Software Ver 4.0.2 (R Core Team, 2020).

Results

The SIN COVID Survey Database [1] included data collected from 365 Italian nephrology centers and involved 30,821 HD, 4139 PD, and 25,481 kidney TX patients, for a total of 60,441 RRT patients.

On April 23rd, 2020, the cumulative incidence of SARS-CoV-2 positive cases was 3.55% (95% CI 3.34–3.76) in HD patients, 1.38% (95% CI 1.04–1.78) in PD patients and 0.86% (95% CI 0.75–0.98) in TX patients [1]. The distribution of SARS-CoV-2 cases among RRT patients was positively skewed (skewness = 5.53), with an average of 3.7 positive patients per center but with a variance of 88.7 (min 0; max 105).

Multilevel modeling of SARS-CoV-2 positive cases

Multilevel negative binomial regression predicted SARS-CoV-2 positivity better than multilevel Poisson and zero-inflated regression (Table 1). The inclusion of random effects at the second level significantly improved the model fit (χ2 = 3.30, p = 0.034). Among candidate second level variables, the model with provinces as cluster-specific random effect had the lowest AIC (Table 1).

Table 1

Model selection for prediction of SARS-CoV-2 infection among HD patients

Model	Second level	Third level	AIC	BIC
Poisson regression			1540.02	1571.16
Zero-inflated Poisson regression			1475.54	1510.57
Mixed-effects Poisson regression	Province		1221.71	1256.74
Negative binomial regression			1108.27	1143.29
Mixed-effects negative binomial regression	Region		1104.78	1143.70
Mixed-effects negative binomial regression	Province	Region	1102.63	1145.44
Mixed-effects negative binomial regression	Province		1100.12	1139.03

AIC Akaike information criterion, BIC Bayesian information criterion

The independent factors (Table 2) associated with increased rates of SARS-CoV-2 positive HD patients were: the geographical latitude of the center (IRR 1.34, 95% CI 1.20–1.51), the positivity rate in the contextual general population (IRR 1.002; 95% CI 1.001–1.003), the count of SARS-CoV-2 positive healthcare workers (IRR 1.09; 95% CI 1.02–1.17), and the test-all policy for both patients and healthcare workers as an interaction term (IRR 5.94, 95% CI 3.36–10.45). Conversely, the number of days since the lockdown date had a protective effect (IRR 0.95, 95% CI 0.94–0.98). The second level effect was significant, and the contextual effect explained, on average, 10% of the total variability [VPC = 0.10 (range 0.002–0.210)]. The net effect of the province as second level variable, holding fixed the first level factors, is shown in Fig. 1 and mapped in Fig. 2.

Table 2

Independent factors associated with the rate of SARS-CoV-2 positive cases in hemodialysis patients

Factors	Units	IRR	Z	p value	95% CI
Geographical latitude	1°	1.34	6.28	< 0.0001	1.20–1.51
Number of positive cases in the dialysis centers’ province	100 cases	1.002	2,60	0.009	1.001–1.003
Positive HCWs	1 case	1.09	3.41	0.001	1.02–1.17
Testing all HCWs	Yes vs No	1.37	1.30	0.194	0.85 -2.31
Testing all patients	Yes vs No	1.67	1.04	0.297	0.64–4.38
Testing all HCWs and patients	Yes vs No	5.94	6.12	< 0.0001	3.36–10.45
Days in lockdown	1 day	0.95	− 5.16	< 0.0001	0.94–0.98
Second level random effect
Variance: 0.283 (SE:0.19)		VPC = 0.104		Range 0.02–0.21

IRR incidence rate ratio, CI confidence interval, VPC variance partition coefficient, HCWs healthcare workers

Fig. 1

Predicted post-estimation counts with 95% confidence intervals by province and ranked from the lowest to the highest. Although the variability is high, only three provinces are significantly different from zero

Fig. 2

Predicted post-estimation count due to the general contextual effect at the province level. A negative value indicates that the contextual effect reduces the expected count, a near null value indicates that the contextual effect is negligible, a positive value suggests an increase in the expected count due to the contextual effect

The SARS-CoV-2 positivity rate in PD patients was associated with the cumulative incidence of positive cases in the province general population (IRR 1.10, 95% CI 1.05–1.20), and with geographical latitude of the centers (IRR 1.42, 95% CI 1.18–1.70) (Table 3). We found no formal evidence for higher-level effects, and the final estimates obtained from a conventional negative binomial model are summarized in Table 2. For TX patients we found the same associations, again with no evidence of a significant contextual effect (Table 4).

Table 3

Independent factors associated with the rate of SARS-CoV-2 positive cases in peritoneal dialysis patients

Factors	Units	IRR	Z	p value	95% CI
Number of positive cases in the dialysis centers’ province	100 cases	1.10	3.07	0.002	1.05–1.20
Geographical latitude	1°	1.31	3.03	0.002	1.10–1.56

IRR incidence rate ratio, CI confidence interval

Table 4

Independent factors associated with the rate of SARS-CoV-2 positive cases in transplant patients

Factors	Units	IRR	Z	p value	95% CI
Number of positive cases in the dialysis centers’ province	100 cases	1.30	8.32	< 0.001	1.20–1.33
Geographical latitude	1°	1.42	3.72	< 0.001	1.18–1.70

IRR incidence rate ratio, CI confidence interval

Classification and regression tree

CART analysis substantially confirmed the predictors identified by the negative binomial model in the following order of importance: geographical latitude of the centers, cumulative incidence of SARS-CoV-2 positive cases in the general population, number of test-positive healthcare workers, and days elapsed since lockdown date. The maximally pruned model generated three rules (Fig. 3):

Fig. 3

The decision tree shows the rules and split points to estimate the number of SARS-CoV-2 positive patients. The first row in each box shows the estimated number of cases, the second row the number of patients and the number of centers, the third row the percentage of cases covered

There were fewer than 1 SARS-CoV-2 positive HD patients per center when the cumulative incidence of positive cases in the province was below 229 per 100,000,

There were nearly 5 SARS-CoV-2 positive HD patients per center when the province cumulative incidence of positive cases was greater than or equal to 229 per 100,000 and there were fewer than 3 test-positive healthcare workers in the center,

There were nearly 15 SARS-CoV-2 positive HD patients per center when the cumulative incidence of positive cases in the province was greater than or equal to 229 per 100,000 and there were at least 3 test-positive healthcare workers in the center.

Multilevel modelling of case fatality

Admission to ICUs directly influenced the infection fatality proportion in SARS-CoV-2 positive HD patients (IRR 1.04; 95% CI 1.01–1.07), with no evidence of clustering in provinces or regions (χ2 = 0.30, p = 0.58). In PD and TX patients, fatality proportions in SARS-CoV-2 positive cases were associated with the number of positive cases in the province general population (IRR 1.003; 95% CI 1.002–1.004 for both TX and PD), again with no evidence of clustering in provinces or regions (χ2 = 0.02, p = 0.89. and χ2 = 0.04, p = 0.84, respectively).

Discussion

The SIN COVID-19 survey was conducted during the exponential phase of the pandemic in Italy. At that time, there were about 190,000 positive cases in Italy [7], and among them, 1093 were HD patients, 57 PD patients, and 218 kidney transplanted patients. The survey provided three main findings [1]: an almost tenfold higher proportion of SARS-CoV-2 positive subjects among RRT patients than in the general population (2.26 vs 0.3%) with a significant difference between HD, PD and TX patients (3.55%, 1.38% and 0.86%, respectively); a sizable variability among Italian regions with a north to south decreasing gradient; a higher case fatality rate among RRT patients than in the general population. Beyond these descriptive findings, we undertook an extended analysis cycle to identify the independent factors associated with SARS-CoV-2 spread and subsequent death.

We initially classified the factors associated with a higher likelihood of SARS-CoV-2 infection in RRT patients into three groups. First, factors shared with the general population, such as the center location and geographical latitude, the local cumulative incidence, and the length of lockdown period. Second, factors related to the center, such as the number of SARS-CoV-2 positive healthcare workers and the testing policy adopted. Third, neighborhood-level contextual factors, such as the organizational characteristics of the local health authorities.

Main factors and contextual effect

Through April 23rd, 2020, significant spread occurred following a north to south gradient in Italy, which is in keeping with a recent study that included eighty-eight countries reporting a statistically significant correlation between country latitude and number of COVID-19 cases [13]. The behavior of a seasonal respiratory virus, influenced by temperature, relative humidity, sunlight exposure [14], and air pollution [15] that is typical of the industrialized north of Italy [16-18], may in part explain such a pattern. On the other hand, the national lockdown measures effectively reduced infection spread, and resulted in a considerable drop in air pollution levels [19]. Hence, the latitude gradient may actually describe the effects of the lockdown on the dynamics of the virus spread. Indeed, the first cases were observed in the northern regions, where the infection rapidly spread, but the implementation of national lockdown measures directly affected the further southwards diffusion of the virus.

The number of positive cases observed in each center is associated with SARS-CoV-2 spread within the province general population. Such finding was somewhat expected because it reflects the prior probability that RRT patients will be infected before any other specific characteristics or condition is considered [20].

The days in lockdown acted as a protective factor and showed that quarantine measures effectively reduce the number of positive cases and infection spread in vulnerable patients. Interestingly, CART analysis identified a split rule for lockdown duration over 42 days, corresponding to three consecutive 14-day quarantines. This finding is indeed interesting because it suggests that shorter lockdown periods may not yield the expected benefit.

Factors related to the center

The greater the number of positive healthcare workers in a center, the higher the likelihood of there being SARS-CoV-2 positive patients. This does not necessarily mean that healthcare workers directly infect patients, but it is rather a proxy of the virus spread within the center. Notably, through CART analysis we found that a threshold of three infected healthcare workers can distinguish centers at low risk from those at high risk of infection spread when the cumulative incidence in the general population exceeds 230 cases per 100,000 inhabitants.

The test-all policy is associated with increased identification of SARS-CoV-2 positive cases, especially when it is adopted for both healthcare workers and patients. An unrestricted testing policy can effectively identify positive patients and contain the contagion, especially if the population positivity rate is on the rise in the surrounding territory [21] because it allows to identify asymptomatic subjects. It is however worth pointing out that the testing policy was modified in a large number of centers during the course of the epidemic. Initially, it was limited to symptomatic cases in most centers, but in keeping with the evolving literature and recommendations, tests were progressively extended to all patients and healthcare workers [22, 23]. Taken together, our findings strongly support the universal testing policy, advocated by countries that first faced the epidemic [24], which is still the topic of a lively debate in most European countries.

Contextual factors

Besides the independent effect of the fixed component of the model, our analysis provides evidence of a general contextual effect. On average, about 10% of the variation in the number of positive cases among centers is due to contextual differences. Hence, the inherent characteristics of geographical and administrative areas to which the centers belong (the context) explain a sizable part of the variation in positivity rates. The areas administered by Italian Local Health Authorities largely coincide with the province territory. Therefore, center policies and resource allocation within the same province is determined by the same Authority, which may explain the observed discrepancies among strategies and timeliness of interventions adopted throughout the country. For instance, in some provinces, lockdown dates, case isolation and management, and directives for containment measures and testing varied according to locally determined policies. Centralization of HD patients, common in China [25], was generally discouraged in Italy, but in some provinces, it was promoted and thus most likely influenced the center-specific rates.

We did not find evidence for a contextual effect among PD and TX patients and they seem to share the same risk factors of the population in the local area. Therefore, promoting home dialysis and kidney transplantation, when possible, is a valuable and effective containment measure in an epidemic context. In this setting, the contextual effect especially explores the set of the organizational characteristics that act beyond the fixed predictors of the outcome, but it does not identify each specific characteristic. Further analyses should be performed to establish the best organizational set-up.

Differences among RRT modalities

Unlike in-center HD patients, the proportion of SARS-CoV-2-positive PD and TX patients depends solely on the density of cases in the general population and the geographical latitude. It is highly conceivable that patients on home-based therapy, with a reduced schedule of center access and a higher likelihood of adherence to general prevention measures, might be less exposed to infection. In other words, these patients share the same degree of risk as the general population in each local area with no evidence for contextual effects.

Unfortunately, our survey data do not allow to investigate the increased immunosuppression-related vulnerability of TX patients. As stated in the survey questionnaire [1], we collected only information concerning the number of total and positive TX patients at the center level and not their individual data, such as transplant history, blood chemistry and immunosuppressive therapy. Nonetheless, it can be argued that any person, be they immunocompetent or immunosuppressed, who comes into contact with a formerly unknown virus has the same susceptibility to the infection, but the consequences are likely to be more severe in immunosuppressed organ transplant recipients [26]. The preliminary evidence reviewed by Thng [27] supports this hypothesis, although further prospective studies are needed to provide conclusive data.

Although RRT patients share the same risk of infection as the general population, only HD patients show additional risks associated with frequent contact with other possibly infected patients and healthcare workers. Clearly, our findings do not provide direct evidence of SARS-CoV-2 spread within the centers, but they do support this hypothesis. Unfortunately, a direct comparison with other groups of patients, regularly treated in hospital facilities for chronic diseases, is hardly feasible. No other patient population undergoes life-saving treatments delivered in the hospital for at least 12 h a week, in the same room, with the same fellow patients for months or years. Perhaps patients undergoing chemotherapy, radiation or physical therapy may somewhat resemble HD patients, but most of these treatments were either re-scheduled or halted during the critical phase of the pandemic [28, 29], consequently they cannot be considered a concurrent “control” group.

RRT patients and the general population

The higher proportion of SARS-CoV-2 infection in RRT patients as compared to the general population is an arguable issue. It is possible that the cumulative number of positive cases in the general population provided by national official statistics may have underestimated the actual case density at the time of the survey. As a matter of fact, the study by Lavezzo et al. [21] shows that the cumulative incidence of SARS-CoV-2 infection in the population in one of the first villages in which the spread occurred in Italy was 2.6% and the preliminary data of a survey on SARS-CoV-2 seroprevalence in Italy [30] show that it was 2.5%. These findings seem numerically consistent with the cumulative incidence observed in RRT patients, but they must be interpreted cautiously because biological samples were collected in a small municipality in the north-eastern part of the country [21] and the seroprevalence was estimated to be less than 1% of the Italian population [30].

SARS-CoV-2-positive RRT patients show higher fatality rates than the general population [31], regardless of treatment and latitude. In HD patients, there is a link with the severity of the disease, expressed by number of ICU admissions, while for PD and TX the fatality rate is associated only with the positivity rate in the general population.

Action thresholds

Through CART analysis we identified a simple and interpretable set of rules to support decision making. The decision tree is a well validated supervised machine learning algorithm used for the classification and regression application. The pathway illustrated in the decision tree (Fig. 3) combines the major determinants of the within-center SARS-CoV-2 positivity rate. Namely, the incidence in the general population, the area of residence, and a proxy of virus spread inside the center (i.e. number of infected healthcare workers). Collectively, the decision rules offer a simple data-driven response tool to promptly act when the risk of infection for the patients becomes critical. Furthermore, our results can help policymakers and clinicians’ decisions-making during the current and future pandemics. Decisions are typically made in an evolving scenario and modeling becomes a valuable tool to provide guidance for establishing prevention strategies and infection control.

Strength and limitations

The main strength of the present survey is that it covers the majority of nephrology centers in a nation that was nearly overwhelmed by the COVID-19 pandemic, but it must be evaluated in the context of its limitations. The findings should be interpreted cautiously since the study database was based on aggregate data, thus making predictive inference somewhat imperfect, despite being based on large numbers. Furthermore, we were not able to stratify patients and outcomes by sex, age, and/or comorbidities because such data were not requested in the survey questionnaire.

Conclusions

To the best of our knowledge, this is the largest nationwide study on COVID-19 in RRT patients carried out to date [1, 22, 23], and the present analysis provides answers to the questions posed in the Introduction. The local number of SARS-CoV-2 positive cases in the general population is likely the most important determinant of the infected proportion among RRT patients, but for HD patients other factors strictly related to in-center models of care contribute to explain the higher frequency of positive cases. Namely, the number of infected healthcare workers, as a proxy of the virus spread within the dialysis center, and the testing policy adopted by the center, especially when it is broadened to include all patients and healthcare workers. These factors explain most of the heterogeneity observed in the descriptive study [1], while about 10% of the remaining variability can be attributed to the contextual organization of local healthcare programs. Finally, our study provides a novel set of decision rules based on thresholds derived from easily accessible and continuously monitored data such as the infected proportions in the local general population and within the dialysis center. If and when action thresholds are crossed, clinicians and policymakers may use them as an indication to quickly adopt more stringent containment measures beyond those suggested by current guidelines [31-33]. Future surveys should use a multi-wave longitudinal design and include more detailed information on both center-level and patient-level variables to explain a greater proportion of the variability.