A dual-center, retrospective, cohort study was performed at the departments of Dermatology of the University Medical Center Groningen (UMCG) and the Erasmus University Medical Center Rotterdam (EMC) in the Netherlands. Inclusion criteria were patients ≥ 18 years with a physician-verified diagnosis of HS and treatment with oral retinoids for HS between 2011 and June 2021. Patients with Down’s syndrome (n = 4) or patients who were treated with retinoids for a different indication (n = 10) were excluded. Because of the retrospective nature of this study, exemption from reviewing was granted by the medical ethical board from the UMCG.

The patient data were extracted from electronic patient files: sex, age at onset of HS symptoms, age at start of retinoid treatment, body mass index (BMI) at start of treatment, smoking status at start of treatment, comorbidities, previous or active acne vulgaris or pilonidal cyst, family history of HS, presence of widespread comedones, phenotype, and the treating center. Phenotypes were categorized as previously published, and stratified accordingly to the regular type, frictional furuncle type, scarring folliculitis type, and conglobata type [9, 15]. Both phenotype and presence of widespread comedones were assessed based on photographic images by three physicians trained in the field of HS (KB, PA, KD) who were blinded to treatment. In the event of discrepancies, the patient in question was discussed among the physicians and a consensus was formed. If no sufficient photographic material was present to reach a conclusion, as was the case in a handful of patients, the phenotype and presence of widespread comedones were extracted from the patient charts. Disease severity was assessed using the International Hidradenitis Suppurativa Severity Score System (IHS4), and the (refined) Hurley score based on patient charts [16,17,18].

Type of therapy (acitretin or isotretinoin), treatment regimens, start and stop dates, reason for discontinuation, previous treatment with retinoids, concomitant medication (topical or systemic), and side effects were listed. Surgical procedures for HS during treatment were recorded and included deroofing, skin-tissue sparing excision with electrosurgical peeling (STEEP), and limited or wide excision surgery.

Normality was assessed using the Kolmogorov-Smirnov test, showing all data to be non-normally distributed. Continuous variables were presented as median and interquartile range (IQR). Categorical data were reported as the number of patients and percentage (n, %). Mann–Whitney U tests were used to analyze continuous data. For categorical variables, the Chi-squared test was used.

Multiple imputation (MI) treats the observed varying percentage of missingness across participants, using the mice package version 3.13.0 on R version 1.4.1103. The data set was tested for the assumption of missing completely at random (MCAR; Little’s MCAR test p = 0.34 > 0.05). BMI, smoking status, ethnicity, (refined) Hurley stage, IHS4 categories, family history, disease duration, and age of onset were imputed through fully conditionally specified models. The continuous IHS4 was deemed not suitable for imputation and was removed from the analysis. Random Forest (for categorical diagnostic outcomes) and Classification and Regression Tree (the rest of the variables for imputation) methods were used to construct the imputation model. The analytical results were pooled based on the Rubin’s rules using in-house developed script (Naderi, et al. under review). For details on the imputed data see Supplementary Methods (appendix S1 in the electronic supplementary material [ESM]).

Per retinoid, the general drug survival was analyzed using Kaplan–Meier survival curves and compared using Cox regression analyses. Separate drug survival curves were generated for discontinuation due to ineffectiveness, side effects, and remission. In the analyses, patients were censored at the date of their last visit, when lost to follow-up, or when the retinoid was stopped for reasons other than ineffectiveness, side effects, or remission. If a patient discontinued treatment due to a combination of ineffectiveness and side effects, it was taken into account as having an event in the general drug survival analyses, ineffectiveness analyses, and side effect analyses. Moreover, selected baseline variables were assessed as possible predictors for longer drug survival time in univariate Cox regression analysis: sex, age, BMI, treating center, disease duration, smoking status, previous retinoid use, widespread comedones, previous or active pilonidal cysts, previous or active acne vulgaris, categorical IHS4, (refined) Hurley stage, concomitant comedication, surgery during treatment, and phenotype. To assess the effect of concomitant immunosuppressive medication and surgery during treatment on the identification of possible predictors for longer drug survival, post-hoc Cox regression sub-analyses were performed excluding patients with concomitant immunosuppressive medication and surgery during treatment. Although classically Cox regression analyses define hazard ratio (HR) as a risk estimate per unit time to reach the outcome, here we used it to identify predictors for a longer treatment survival. Assumptions of the Cox regression analyses were tested using the time-dependent explanatory variables; if a p < 0.05 was found, the assumptions were not met and left out of the analyses. In the statistical analyses, differences with two-sided p values < 0.05 were regarded as statistically significant. Statistical analyses were performed using IBM SPSS Statistics for Windows, version 25 (IBM Corp., Armonk, NY; USA).