This study was divided into three parts. The first part was the preparatory process before using the administrative database (DB); we validated ERT patients using an identification method and validated a risk-adjusted mortality model incorporating important clinical variables for ischemic stroke (Figure 1A and B). The second part was to examine VORs by analyzing the association between hospital ERT volume and risk-adjusted outcomes derived from summing the individual-level expected outcomes using the administrative DB (Figure 1C). The third part was to calculate the hospital volume threshold (Figure 1D).

Three clinical outcome variables, namely 30-day mortality, 30-day readmission, and postprocedural intracranial hemorrhage (ICH), were selected based on literature review [21-23]. After identifying other variables for building risk adjustment models for three outcomes from the administrative DB for hospital profiling [22], we requested and purchased these variables and other related claims data from the Health Insurance Review and Assessment Service (HIRA) in Korea. Next, we verified the administrative DB with the real-world data (clinical DB I (Figure 1A) from two university hospitals: Dong-A University Hospital (DAUH) and Pusan National University Yangsan Hospital (PNUYH) [24]. We also considered proxy values to replace the National Institutes of Health Stroke Scale (NIHSS), which is difficult to obtain directly from the administrative DB. We used the stroke severity index (SSI) (Supplementary Table 1) [25]. To achieve this, we obtained the clinical DB II (Figure 1B) and the local administrative claims DB of all ischemic stroke patients from the aforementioned university hospitals. The derivation and validation of the proxy values were performed separately using different DBs.

From the administrative DB in the HIRA (the HIRA DB), risk-adjusted outcome models for each of the three outcomes were derived by using multiple regression analysis [22,26]. Many variables, including modified Charlson comorbidity index (mCCI), and SSI were also analyzed in the model derivation (Supplementary Table 2). The mCCI is one of the variables often used for risk adjustment purposes in the studies using administrative DBs [27]. We measured C-statistics to measure the discrimination power for each of the three models. Next, these models were applied to the HIRA DB, and the outcomes were measured. The relationship between the procedural volume per hospital and the outcomes were then plotted and analyzed to identify VOR (Figure 1C). Based on the results of the VOR analysis, the volume threshold was measured. Finally, to investigate the modified structures and processes, a stratified analysis with hospital characteristics was performed (Figure 1D).

All analyses were conducted with the SAS Enterprise Guide statistical software package version 7.13 (SAS Institute Inc., Cary, NC, USA).

To gain access to HIRA DB, we submitted an Institutional Review Board approval document from PNUYH to HIRA. There were many datasets including diagnostic codes of cerebrovascular diseases (such as cerebral ischemic stroke and subarachnoid hemorrhage) and various procedural codes (for instance, mechanical and chemical thrombectomy, Levin tube insertion, and coil embolization of intracranial aneurysm.) From these claims data, 12,013 patients with cerebral AIS who underwent ERT from July 2011 to June 2016 in 162 hospitals were selected. Information about patient and hospital were accessed in a deidentified state.

Clinical DB I was a prospective clinical cohort for patients with AIS undergoing ERT in DAUH and PNUYH. The enrollment period was the same as in the HIRA DB. Clinical DB II was a cohort of all patients with AIS in addition to ERT patients during the same enrollment period and linked with the local claims database of the same patients (Figure 1A and B).

Although the HIRA DB was used in a de-identified state, it did provide geographical information (city and district). Combining this information with the relative volume of ERT, we produced the DAUH and the PNUYH cohorts from the HIRA DB. To validate these cohorts, a 2×2 contingency table was created, and Clinical DB I was defined as the gold standard. Pearson’s chisquare test and Student's t-test were used to compare the administrative and clinical databases (Methods section in Supplementary appendix [24,25,27-31]).

To derive the SSI, the association between NIHSSs of patients with AIS and SSI were analyzed for all patients with AIS of the DAUH (prediction model). Subsequently, the model was validated with the PNUYH data. The correlation coefficients were 0.70 for the prediction model and 0.68 for the validation model (Supplementary Table 3). The correlation coefficients of both models were similar to those in a previous paper [25] and usually considered relatively high. However, “intensive care unit stay” showed no significant statistical power (Supplementary Table 4). Multiple linear regression analysis was used for creating the SSI model, and Pearson’s chi-square test was used for validating the model (Figure 1B, Methods section in Supplementary appendix).

Baseline characteristics are shown in Table 1. Of the 111 hospitals included in the VOR analysis, the ERT volumes corresponding to half of all patients were performed in 21 hospitals. Age and the proportion of comorbidities appeared to be higher in the first quartile group than in the other groups.

A total of 51 hospitals (31.4%) were excluded in the VOR analysis due to very low hospital volumes. VOR analyses were conducted on the total number of cases for the remaining 111 hospitals, which was 11,745 (97.7% of the original 12,013 cases). The sample was divided into four groups to have roughly the same number of patients in each group. About 2,900 cases from the most voluminous hospitals were categorized as the first quartile (more than 61 cases per year), the second quartile (31 to 60 cases per year), the third quartile (19 to 30 cases per year), and the fourth quartile (less than 18 cases per year) (Table 1).

Over 5 years, 12,013 ERT procedures were performed (Table 1). Some patients (0.009%) had undergone the procedure twice (n=111) or three times (n=2). Of the 162 hospitals, 21 did not perform the ERT procedure between July 2015 and June 2016. The mean volume was 14.8 cases per year, and 59 hospitals (38%) had a volume that was higher than the mean. The overall 30-day mortality rate was 11.6%, 30-day readmission rate was 4.6%, and postprocedural ICH rate was 8.6%.

The results of validating ERT patient identification were a sensitivity of 83.2%, a specificity of 98.8%, and a positive predictive value of 80.2% (Figure 1A). To evaluate the discrepancies, in-depth analyses were performed for the 115 true positives and 42 false negatives of the PNUYH cases (Supplementary Tables 5 and 6).

To solve the multicollinearity problem, other variables included in the SSI were excluded from the univariate analysis, and the statistical model was applied. In the case of mortality, the variables applied to the model were age, SSI, mCCI, sex, hypertension, and atrial fibrillation. The C-statistics of the final mortality rate model was 0.8146, which was higher than the initial target value of 0.8. Similarly, the same modeling was done for the readmission rate, and included variables were SSI, mCCI, age, transient ischemic attack, prior stroke, and female sex. The C-statistic was 0.6052, which was relatively low. For postprocedural ICH rate, individual-level variables, such as age, SSI, intracranial arterial stenosis, myocardial infarction, hypertension, and transient ischemic attack were selected. ‘Teaching status: resident’ was included at the hospital-level. The final model C-statistic was 0.6841 (Supplementary Figures 1-3).

For mortality rate comparison between quartile groups, the odds ratio (OR) showed differences between the 1st and 3rd quartiles, and between the 1st and 4th quartiles in the crude outcome model (P<0.05). However, in the adjustment model, the OR showed differences among the groups from the 2nd to the 4th quartile (1st vs. 3rd, P=0.003; 1st vs. 4th, P<0.001). Similar results were obtained on comparing readmission and postprocedural ICH rates among quartile groups (Table 2). In the distribution analysis, three outcomes of both crude and risk-adjusted models showed negative correlations, but there were no statistically significant results (Supplementary Figures 4-6).

In the between quartile comparison of hospital volume, there was a significant difference in the OR between the 1st and 2nd quartiles for the risk-adjusted mortality rate, and between the 1st and 4th quartile for the risk-adjusted postprocedural ICH rate. Based on this, we concluded that there was a VOR for ERT in patients with AIS, and subsequently estimated the volume threshold.

It was assumed that the threshold would exist between the volume of the boundary point between the 2nd quartile (less than 60 cases per year) and the boundary point of the 4th quartile (more than 19 cases per year). Between 60 and 19 cases per year, we assessed volume threshold using a modified method presented in a previous paper [33]. In this process, a reversed J-shaped curve was created, and the threshold was set to the point at which the odds for mortality approached “1.0” as the total volume falls. Thus, at the highest value of the J-curve, 22 cases per year corresponded to the inflection point, and 24 cases per year corresponded to where the odds fell below 1. Based on this, we set 24 cases per year as the VOR threshold of ERT for AIS (Figure 2).

To identify modifiable ‘structures and processes,’ additional analyses were conducted regarding the two hypothesized hospital-related variables, stroke unit [28], and teaching status [29]. Of the 111 hospitals included in the analysis, 50 hospitals (45%) had stroke units, and teaching status hospitals accounted for 79% of all hospitals. The proportions of hospitals with stroke units and teaching status decreased from the 1st to 4th quartile group, and was statistically significant in the fourth quartile group. Due to dichotomization of the classification according to the value determined by the volume threshold, significant differences were found in the proportion of both stroke unit and teaching hospitals between the two quartiles (Table 3).