This study included case-control samples from the MASTOS study, a population-based case-control study of breast cancer in Cyprus (Loizidou et al., 2008; Hadjisavvas et al., 2010). The study was approved by the National Bioethics Committee of Cyprus, and all participants provided signed informed consent (Approval No, ΕΕΒΚ/ΕΠ/2016/38) (Loizidou et al., 2008). The study includes 990 breast cancer cases and 1,094 age-matched healthy controls. All study participants were over 18 years of age. The average age of diagnosis for cases was 51.5 ± 9.3 standard deviation (sd) years, while the average age at the time of study enrolment for age-matched controls, was 55.7 ± 6.9 sd. Positive family history of breast cancer was reported for 16.2% (n = 158/973) and 8.0% (n = 87/1,091) for cases and controls, respectively. The majority of tumors were invasive carcinomas (92.5%, n = 656/709), estrogen receptor (ER)-positive (75.8%, n = 476/628), progesterone receptor (PR)-positive (61.6%, n = 381/619) and human epidermal growth factor receptor 2 (HER2)-negative (83.9%, n = 490/584). Out of the 990 tumors, 576 (58.2%) had available data for all ER, PR and HER2, of which 11.5% (n = 66/576) had a triple-negative phenotype, 76.6% (n = 441/576) had an ER + phenotype, 62.8% (n = 362/576) had a PR + phenotype and 16.1% (n = 93/576) had a HER2+ phenotype. All study sample characteristics are summarized in Table 1.

We analyzed targeted panel sequencing data on 35 actionable and suspected breast cancer susceptibility genes (Supplementary Table S1) using the BRIDGES panel (Dorling et al., 2021). Results are presented for variants in 34 genes since PPM1D was shown to be associated with breast or ovarian cancer risk, but in low allelic fractions (“somatic mosaicism”) (Ruark et al., 2013). These are potentially due to treatment, thus excluded from further analysis. Details on library preparation, sequencing and bioinformatics analysis including variant calling and quality control were previously documented (Dorling et al., 2021). In-frame insertions/deletions, intronic variants and variants in untranslated regions (UTRs) were not considered in the analysis.

Protein-truncating variants (PTVs) were defined as frameshift insertions/deletions, splicing variants (±2 positions) and nonsense variants as annotated by ANNOVAR (Wang et al., 2010). Splice variants affecting the penultimate exon (except variants in ATM, BARD1, BRCA1, RAD51C, RAD51C and PALB2 for which the truncated protein might still be pathogenic, irrespective of exon skipping; evidence previously documented (Dorling et al., 2021)), as well as PTVs in the last exon of each gene, were excluded from the analysis (Dorling et al., 2021). In addition, six canonical splice variants in BRCA1 (c.594-2A>C, c.4096 + 1G>A, c.4096 + 2T>C, c.4186-2A>G, c.4358-1G>C and c.4358-2del) were excluded, since they are of uncertain clinical significance according to ENIGMA (Evidence-based Network for the Interpretation of Germline Mutant Alleles) classification schemes (Spurdle et al., 2012).

Classifications for rare missense variants (allele frequency of less than 0.001 in gnomAD v2.1.1 non-Finnish European exome samples) were retrieved from ClinVar (https://www.ncbi.nlm.nih.gov/clinvar/, last accessed 11/12/2022), or the ENIGMA BRCA1 and BRCA2 expert panel (https://enigmaconsortium.org/) (Spurdle et al., 2012), along with missense TP53 variant classification based on ACMG/AMP guidelines (Richards et al., 2015), augmented with classifications made using a quantitative classification model which utilizes bioinformatics prediction tools alongside germline to somatic frequency ratios (Fortuno et al., 2019; Fortuno et al., 2021). Rare missense variants were further annotated for functional protein domain location defined by the UniProt database (https://www.uniprot.org/, release 2022_04). Rare missense variants classified as (likely) benign were not considered in the analysis.

Prevalence of protein-truncating and rare missense variants in each gene was tabulated for breast cancer patients and controls. Odds ratios (ORs) with 95% CIs were estimated using Firth’s bias-reduced penalized-likelihood logistic regression for overall, ER-positive and ER-negative breast cancer. For ER-negative cases, we also evaluated the associations for triple-negative breast cancer. Associations were adjusted for age at diagnosis or interview and first-degree family history of breast cancer. Logistic regression p values were estimated using Wald’s test. p values smaller than 0.05 were considered to be statistically significant. Association analysis was performed separately for PTVs and rare missense variants. Associations for rare missense variants were also evaluated according to gene and gene domain for the established breast cancer susceptibility genes (Dorling et al., 2021). Heterozygous and homozygous carriers of variants in a gene were not distinguished as it was not always possible to do so with certainty, and the number of homozygotes was too small for separate analysis.

Cumulative risks of breast cancer in the absence of other events were calculated by combining age-specific odds ratios with population incidence rates for Cyprus (2020) as a baseline, as previously described (Schmidt et al., 2016).

Among 990 breast cancer cases and 1,094 controls, we identified 27 unique PTVs in 4.55% (n = 45/990) of cases and 1.46% (n = 16/1,094) of controls (Table 2, Supplementary Table S2). Of these, 3.54% of cases (n = 35/990) and 0.37% of controls (n = 4/1,094), harbored PTVs in the established breast cancer susceptibility genes ATM, BARD1, BRCA1, BRCA2, CHEK2, PALB2, RAD51C, RAD51D and TP53 (Dorling et al., 2021).

Among the case patients, the highest prevalence of PTVs was observed for BRCA2 (n = 18/990, 1.82%), ATM (n = 9/990, 0.91%), PALB2 (n = 5/990, 0.51%), RAD50 (n = 4/990, 0.40%), BRCA1 (n = 3/990, 0.3%) and FANCC (n = 3/990, 0.3%). We also identified one PTV in each of RECQL, ABRAXAS1 and MUTYH. Among controls, the highest prevalence was observed for RECQL (n = 3/1,094, 0.27%) and FANCC (n = 2/1,094, 0.18%). We also detected one PTV in each of ATM, BRCA1, BRCA2, CDH1, CHEK2, FANCM, GEN1, MRE11, MUTYH, RAD50 and RINT1 (Table 2). The founder BRCA2 c.8756delG PTV was detected in 1.01% of breast cancer cases (n = 10/990) and 0.09% of controls (n = 1/1,094) which corresponds to more than half of the BRCA2 PTVs observed among cases (55.56%, n = 10/18) and the only BRCA2 PTV identified among controls.

We report evidence of PTV association with overall breast cancer risk for BRCA2 and ATM with adjusted ORs of 9.75 (1.81-52.69 95% CI, p-value = 0.0081) and 7.61 (1.32-43.97 95% CI, p-value = 0.023) (Figure 1; Table 2, Supplementary Table S3). Association with overall breast cancer is also depicted for PTVs in PALB2 with adjusted ORs of 16.14 (0.85-305.14 95% CI, p-value = 0.064), although not statistically significant (Figure 1; Table 2, Supplementary Table S3). For other genes included in the panel, the evidence for an association between overall breast cancer and PTVs did not reach statistical significance (Figure 1; Table 2, Supplementary Table S3).

Breast cancer cases with PTVs in BRCA1 and BRCA2 had a relatively younger age at diagnosis (46.7 ± 10.4 sd) compared to cases without PTVs in BRCA1 and BRCA2 (51.6 ± 9.3 sd, Wilcoxon rank sum test p-value = 0.046) or cases with PTVs in genes other than BRCA1 and BRCA2 (53 ± 10.9 sd, two-sample t-test p-value = 0.051). The associations with breast cancer overall for both participants diagnosed or recruited at an age earlier or later than 50 years did not reach statistical significance, except for BRCA2 PTV-carriers who were at an increased risk of the disease if diagnosed at an age later than 50 years (Supplementary Table S4). Analysis by family history of breast cancer (first-degree relatives) revealed that the prevalence of PTVs among patients with a family history of breast cancer (8.23%, n = 13/158) was about two times the frequency reported among patients with no family history of breast cancer (3.56%, n = 29/815; p-value = 0.015). To investigate the influence of a family history of breast cancer, analysis was conducted separately for cases with and without a first-degree relative with breast cancer. Among participants with family history of breast cancer, PTVs in ATM (p-value = 0.002), BRCA2 (p-value = 0.005), BRCA1 (p-value = 0.046), PALB2 (p-value = 0.039) and RAD50 (p-value = 0.044) were significantly associated with increased risk of breast cancer. The exclusion of breast cancer patients with a family history of breast cancer had a minor impact on PTV association with overall breast cancer risk (i.e., BRCA2 adjusted OR of 8.74; 1.57-48.75 95% CI, p-value = 0.013) (Supplementary Table S5).

Protein-truncating variants in BRCA2 were strongly associated with an increased risk of ER-positive and triple-negative breast cancer with adjusted ORs of 10.18 (1.77-58.55 95% CI, p-value = 0.009) and 10.14 (1.08-94.84 95% CI, p-value = 0.042), respectively. Protein-truncating variants in BRCA1 were strongly associated with an increased risk of ER-negative breast cancer with adjusted ORs of 20.73 (2.04-210.69 95% CI, p-value = 0.0104). Protein-truncating variants in ATM were associated with an increased risk of both ER-negative and ER-positive disease with adjusted ORs of 12.55 (1.33-118.68 95% CI, p-value = 0.027) and 7.47 (1.09-51.23 95% CI, p-value = 0.041), respectively. Finally, PTVs in the PALB2 gene were associated with an increased risk of ER-negative breast cancer with adjusted ORs of 44.23 (1.78-1098.24 95% CI, p-value = 0.021) (Supplementary Table S3).

Odds ratios decreased significantly with increasing age for BRCA2 with ORs of 0.94 (0.89-0.99 95% CI, p-value = 0.013) (Supplementary Table S6). Estimated cumulative risks of breast cancer in the absence of other events were calculated by combining age-specific odds ratios with population incidence rates for Cyprus (2020) (Figure 2). For carriers of PTVs in BRCA2 and ATM the estimated cumulative risks by 80 years of age exceeded the 30% threshold for high risk, as defined by the NICE (National Institute for Health and Care Excellence) surveillance screening guidelines (NICE, 2019).

Among 990 breast cancer cases and 1,094 controls, we identified 380 unique rare missense variants in 50.40% of cases (n = 499/990) and 41.77% of controls (n = 457/1,094) (Supplementary Table S2). Of these, 147 unique rare missense variants were detected in the established breast cancer susceptibility genes, among 195 out of 990 cases (19.70%) and 188 out of 1,094 controls (17.18%). None of the detected rare missense variants were classified as pathogenic.

Rare missense variants in BRCA1 were associated with an increased risk of overall breast cancer with an adjusted OR of 2.4 (1.14-5.08 95% CI, p-value = 0.022) (Figure 3). In the subtype-stratified analyses, rare missense variants in BRCA1 were associated with an increased risk of both ER-negative and ER-positive disease with adjusted ORs of 3.56 (1.15-11.01 95% CI, p-value = 0.028) and 2.77 (1.21-6.38 95% CI, p-value = 0.016), respectively (Supplementary Table S7, Figure 3). Increased risk of triple-negative breast cancer was observed for carriers of rare missense variants in PALB2 with an adjusted OR of 5.79 (1.32-25.29 95% CI, p-value = 0.020, Figure 3).

For missense variants in aggregate, the evidence for an association between overall breast cancer and any protein domain did not reach statistical significance. However, there was a significantly increased risk of ER-negative breast cancer associated with variants located within the BRCA1 Zinc-Finger, ATM FAT, PALB2 WD1 and CHEK2 FHA domains. In addition, rare missense variants within the ATM PI3L/PI4K catalytic domain demonstrated an association with an increased risk of ER-positive disease. Furthermore, rare missense variants within PALB2 domains WD1, domain for stimulation of POLH DNA synthesis and domain for interaction with RAD51, BRCA2 and POLH were significantly associated with increased risk of triple-negative breast cancer. Detailed results from the domain-specific analysis are shown in Supplementary Table S8.