We analyzed 54 IBC samples and 102 non-IBC TCGA samples (Additional file 1: Table S1). All cases were ductal type and from women with AJCC stage III-IV à. In the IBC group, the median patients’ age at diagnosis was 48 years (24–79) and 64% were non-menopausal; the pathological grade was 3 in 70% of cases and the molecular subtypes were HER2+ in 41% of cases, HR+/HER2− in 37%, and TN in 22%. As expected, IBCs were associated with younger patients’ age, more frequent non-menopausal status and HER2+ and TN subtypes, than non-IBCs.

WES analysis identified 5576 somatic mutations in 4200 genes in IBCs and 6749 somatic mutations in 4839 genes in non-IBCs (Additional file 2: Table S2, Additional file 3: Table S3). The median number of somatic mutations per sample was not significantly different between both groups (62.5 in IBCs (3–942) vs. 50 in non-IBCs (4–442), p = 0.349, Wilcoxon test). The median TMB was not different between IBCs and non-IBCs (1.24 mutations/MB (0.05–21) vs. 1.32 (0.10–12) respectively, p = 0.294, Wilcoxon test; Additional file 7: Fig. S1A), even after adjustment on the molecular subtypes (p = 0.433, logit function). However, after adjustment for age differences and for molecular subtypes, the TMB was approximately 20% higher in IBC (OR = 1.18; P = 0.05). Four percent of IBC samples presented a high TMB (> 10 mutations/MB) vs. 1% of non-IBC samples (p = 0.275, Fisher’s exact test). The TMB was higher in the TNBC subtype (Additional file 7: Fig. S1B) than in the HR+/HER2− subtype in both IBCs (p = 0.115) and non-IBCs (p = 0.062). The median number of tumor neoantigens per sample was not significantly different between IBCs (47, range 9 to 97) and non-IBCs (34, range 3 to 574) without (p = 0.171, Wilcoxon test) and with (p = 0.963, logit function) adjustment on the molecular subtypes.

Among the somatic mutations, 96% were single nucleotide variants (SNVs) and 4% were insertions/deletions (indels) in IBCs versus 90% and 10% respectively, in non-IBCs (p = 2.21E−42, Fisher’s exact test). The difference remained significant after adjustment on the molecular subtypes (p = 8.41E−39, logit function). The percentage of non-silent mutations was lower in IBCs (74%) than in non-IBCs (78%; p = 8.99E−07, Fisher’s exact test), even after adjustment on the molecular subtypes (p = 4.28E−07, logit function). If we consider the SNVs only, the percentage of non-silent mutations was also lower in IBCs (73%) than in non-IBCs (75%; p = 4.49E−03, Fisher’s exact test), and the difference remained significant after adjustment on the molecular subtypes (p = 3.23E−03, logit function). Among the non-silent SNVs, the percentage of missenses was higher (95% vs. 93%) and the percentage of nonsenses was lower (5% vs. 7%) in IBCs than in non-IBCs (p = 2.74E−04, Fisher’s exact without adjustment for the molecular subtypes and p = 4.13E−04, logit function with adjustment).

A total of 195 out of 5576 mutations (3.5%) were defined as driver mutations (TIER1-TIER2) by Cancer Genome Interpreter (CGI) in IBCs versus 357 (5%) in non-IBCs (Wilcoxon test: p = 4.60E-02 without adjustment for the molecular subtypes, but p = 0.829 after adjustment). They concerned 117 genes in IBCs (Additional file 2: Table S2), including classical driver genes of breast cancer, such as TP53, PIK3CA, MAP2K4, GATA3, or KMT2C. The 31 genes mutated in at least 2 IBCs are shown in Fig. 1. Forty-eight percent of them are included in the 93-gene list of driver genes defined by Nik-Zainal et al. in non-IBCs [20]. The most commonly mutated gene was TP53 (54% of samples), followed by PIK3CA (22%). All 117 genes displayed similar mutation frequency between IBCs and non-IBCs in our series (p > 0.05; Fisher’s exact test).

Distribution of alterations of the top 31 genes mutated in IBCs. Oncoprint of the top 31 genes mutated in at least two IBC samples. Top: Number of smatic mutations in each sample. IHC-based molecular subtypes and IBC/non-IBC groups are color-coded as indicated in the legend. Bottom: somatic gene mutations color-coded according to the legend. The genes are ordered from top to bottom by decreasing percentage of altered IBCs right panel). The percentages of mutation in IBCs and non-IBCs are shown to the right of the Oncoprint

The proportions of base substitutions across SNVs are shown in Fig. 2A, B. In IBCs, the most frequent base change was C > T (average of 48% of substitutions) with respect to single-nucleotide-mutation contexts (Fig. 2A), as observed in non-IBCs (average of 52%). The mutational spectra were similar between IBCs and non-IBCs, except an enrichment in C > G in IBCs (average of 18% vs 14%, p = 2.6E−02, Wilcoxon test; and p = 0.070 after adjustment for the molecular subtypes). The same analysis regarding the tri-nucleotide mutation contexts (Fig. 2B) showed that the most frequent base change in IBCs was G[C > T]G, as observed in non-IBCs; the comparison between IBCs and non-IBCs revealed an enrichment in 11 substitutions and tri-nucleotide contexts in IBCs, notably the G[C > G]G (p = 1.01E−03, Wilcoxon test) and the G[C > T]C (p = 1.42E−02, Wilcoxon test), but none remained significant after FDR correction.

Mutational processes of somatic SNVs in IBCs. A Proportions of base substitutions with respect to single-nucleotide-mutation contexts in IBCs and non-IBCs. B Similar to A but with respect to tri-nucleotide mutation contexts. C Proportions of the most represented COSMIC mutational signatures in the whole population age-related: signature 1; homologous recombination deficiency HRD: signature 3; APOBEC activation: signatures 2 and 13; mismatch repair: signatures 6, 20 and 26; POLE: signature 10). The signatures, IHC-based molecular subtypes and IBC/non-IBC groups are color-coded according to the legend

We then assessed the distribution of the 30 COSMIC mutational signatures. In both IBCs and non-IBCs, the most represented signatures were, as expected, the signatures 1 (age-related), then 2 and 13 (APOBEC activation), then 3 (homologous recombination deficiency), and then 6, 20 and 26 (mismatch repair) (Fig. 2C). In IBC, positive correlations existed between the abundances of signature 2 (Additional file 8: Fig. S2A) and of signature 13 (Additional file 8: Fig. S2B) and higher TMB (p = 4.0E−04 and p = 5.6E−03 respectively, Wilcoxon test), and between the abundances of signature 3 (Additional file 8: Fig. S2C) and higher HRD score (p = 1.2E−03, Wilcoxon test), The comparison with non-IBCs identified only one signature differentially represented between both groups: the signature 1 was less frequent in IBCs than in non-IBCs (p = 1.69E−02, Wilcoxon test), even after adjustment on molecular subtypes (p = 3.91E−02, logit function; Fig. 2C).

Given the enrichment in C > G transversions in IBC, we set out to identify possible new mutational signatures. Using non-negative matrix factorization (NMF), six signatures were extracted from our data. Four of these were associated with mutational processes associated with APOBEC activity (N = 2), mismatch repair (N = 1), and homologous recombination deficiency (N = 1). The two remaining signatures (i.e. SBSA and SBSB) have an unknown etiology and were not recovered from the non-IBC mutational profiles. SBSA is characterized by a rather flat profile, whereas SBSB shows dominance of C > T transitions (Additional file 8: Fig. S2D).

Figure 3 shows the frequency plots of low- or high-level CNAs. In non-IBCs, the most frequently gained regions were on 1q, 8q, 11q, 17q and 20q chromosomal arms, whereas the regions frequently lost were on 8p, 11q and 16q. Globally, visual inspection did not reveal obvious differences between IBCs and non-IBCs in terms of altered regions and of frequencies of alterations. GISTIC analysis of IBC samples identified 28 chromosomal cytobands significantly (q < 0.25) gained/amplified (total length of 92 Mb) and 18 chromosomal cytobands significantly (q < 0.25) lost/deleted (total length of 854 Mb). The gained/amplified cytobands comprised 725 genes including 8 defined as driver alterations by CGI: HER2, CCND1, MYC, EGFR, PIK3CA, FGFR2, MDM4, and AKT3, as well as FGFR1, ZNF703. The two most significant gained/amplified cytobands (17q12 and 11q13.3) were regions classically amplified in breast cancer. As expected, all HER2-negative IBC tumors (by IHC/FISH) had no HER2 gain/amplification, whereas 19 out of 22 HER2-positive IBC tumors had HER2 gain/amplification. The three discordant tumors have no HER2 mutation and likely represent false negatives probably due to sampling bias such as normal tissue contamination or tumor heterogeneity. The lost/deleted cytobands comprised 6,427 genes, including 25 identified as driver genes by CGI, such as NF1, TP53, CDKN1A, ATM, STK11, BAP1, ARID1A (Additional file 4: Table S4). We compared the alteration frequencies between IBCs and non-IBCs of genes included in the GISTIC regions gained/amplified in IBCs (Additional file 5: Table S5) and of genes included in the GISTIC regions lost/deleted in IBCs (Additional file 6: Table S6). No gene was more frequently gained/amplified in IBCs than in non-IBCs. Thirty-seven genes (located on 8q21 and including IL7 and HEY1) were more frequently lost/deleted in IBCs than in non-IBCs (p < 0.05; Fisher’s exact test), but none of them remained significant after FDR correction (p > 0.5).

Frequency plots of CNAs in IBCs. Frequencies vertical axis, from 0 to 100%) are plotted as a function of chromosome location for IBCs top) and non-IBCs middle). Vertical lines indicate chromosome boundaries. The CNAs are color-coded as indicated in the legend: gains red), amplifications dark red), losses green), and deletions dark green)

The HRD score was not different between IBCs (median = 28, range 1–99) and non-IBCs (median = 27, range 3–87; p = 0.728, Wilcoxon test; Fig. 4A). By using the classical positivity cut-off (score = 42), 27% of IBC samples were defined as “BRCAness” versus 17% of non-IBC samples, but the difference was not significant (p = 0.195, Fisher’s exact test; Fig. 4B). As expected, this HRD score was higher in the TN subtype than in the HR+/HER2− subtype in both IBCs (p = 0.014) and non-IBCs (p = 3.0E−03, Wilcoxon test; Fig. 4C). Interestingly, it also tended to be higher in TN IBCs than in TN non-IBCs (p = 0.08, Wilcoxon test).

HRD score, heterogeneity index and mutational clonality in IBCs. A Box-plot of HRD score in non-IBC and IBC samples. B Contingency table between HRD score and IBC/non-IBC groups. C Similar to A/, but per molecular subtype. D Box-plot of Heterogeneity H) index in non-IBC and IBC samples. E Contingency table between the tumor heterogeneity status and IBC/non-IBC groups. F Box-plot of the percentages of clonal and subclonal mutations in non-IBC and IBC samples

We also measured the intratumor heterogeneity of tumor samples using SciClone. Heterogeneity index (H-index) was slightly lower in IBCs (median = 0.65, range 0–1.92) than in non-IBCs (median = 0.86, range 0–2.13) samples (p = 0.282, Wilcoxon test; Fig. 4D). Twenty-one percent of IBCs (11/52) vs 33% (32/95) of non-IBCs displayed an H-index superior to 1 (OR = 1.88; p = 0.131, Fisher’s exact test), corresponding to more heterogeneous tumors than the ones with an H-index inferior to 1. Finally, we assessed the percentage of clonal or subclonal non-synonymous mutations in all samples (Fig. 4F). The proportions of clonal mutations were higher in IBCs than in non-IBCs (p = 1.4E−04, Wilcoxon test; p = 3.0E−04 after adjustment upon the molecular subtypes), and consequently subclonal mutations were more frequent in non-IBCs. In IBCs, a positive correlation existed between the abundance of APOBEC signature 13 (Additional file 8: Fig. S2E) and that of subclonal mutations (p = 0.217, Wilcoxon test). When assessing the number of tumor neoantigens in function of the clonality of non-synonymous mutations, we did not observe significant differences in terms of clonal neoantigens (p = 0.980) but the number of subclonal neoantigens tended to be lower in IBC (p = 0.098). In line with this, tumor neoantigens are significantly more often subclonal in non-IBC (p = 0.016), but not in IBC (p = 0.455).

Using the OncoKB database of actionable genetic alterations (AGAs) [42], two levels of clinical evidence (LOE) were distinguished: LOE 1–2 corresponding to standard care therapies, and LOE 3–4 corresponding to investigational therapies. Overall, 72% of IBC samples (39/54) had at least one AGA, versus 68% of non-IBC samples (69/102) when considering all LOE pooled (p = 0.589, Fisher’s exact test; Fig. 5A). Regarding the LOE 1–2 AGAs, 54% of IBCs (29/54) displayed at least one alteration versus 50% of non-IBCs (51/102; p = 0.737, Fisher’s exact test; Fig. 5B). In IBCs, the LOE 1–2 alterations included HER2 amplifications (35%; 19 patients), PIK3CA mutations (22%; 12 patients), and BRCA1 mutation or deletion (6%; 3 patients). For LOE 3–4 AGAs, these figures were 31% in IBCs (17/54) and 36% in non-IBCs (37/102; p = 0.598, Fisher’s exact test; Fig. 5C). Of note, 12/54 IBC samples (22%) had two or more AGAs simultaneously, suggesting potential interest of drug combinations, including two patients with double level 1 PIK3CA mutation. The same was observed in non-IBCs, with 27/102 samples (26%) having two or more AGAs, including 2 patients with double level 1 PIK3CA mutation. The same analysis was done using the ESCAT LOE I–II AGAs. Fifty-five percent of IBCs (30/54) displayed at least one alteration versus 48% of non-IBCs (49/102; p = 0.403, Fisher’s exact test; Fig. 5D). In IBCs, the identified ESCAT LOE I AGAs included HER2 amplification (35%; 19 patients), and PIK3CA mutations (19%; 10 patients). No germline BRCA1/2 mutation, nor MSI status, nor NTRK fusion were identified. The ESCAT LOE II AGAs included PTEN deletion (6%; 3 patients) and AKT1 mutation (2%; 1 patient). No other ESCAT II alteration (ESR1 mutation, HER2 mutation) was identified. In both IBC and non-IBC, profiles of AGAs in PIK3CA and ERBB2 AGAs appeared to be mutually exclusive, although our sample size is too limited to obtain statistical significance.

Percentages of patients with AGAs in IBCs. A Bar-plots of the percentages of patients with IBC and non-IBC displaying at least one OncoKB AGA. The p-value is for the Fisher’s exact test. B Similar to A, but for OncoKB LOE 1–2 AGAs. C Similar to A, but for OncoKB LOE 3–4 AGAs. D Similar to A, but for ESCAT LOE I–II AGAs