HTRA1 expression is associated with immune-cell infiltration and survival in breast cancer

Introduction

Breast cancer (BRCA) is currently the most frequently diagnosed cancer type and one of the leading causes of cancer-related deaths among women worldwide (1). In clinical practice, immunohistochemical staining for estrogen receptor (ER), progesterone receptor (PR), and erb-b2 receptor tyrosine kinase 2 (ERBB2)/human epidermal growth factor receptor 2 (HER2) are utilized as diagnostic markers to categorize BRCA into luminal A, luminal B, HER2-enriched (non-luminal), and basal [triple-negative breast cancer (TNBC)] subtypes (2). Advances in early detection, treatment, and promising therapeutic opportunities have significantly improved the clinical outcomes and survival rates of BRCA patients (3,4). However, the heterogeneity of BRCA subtypes contributes to varying clinicopathological characteristics, treatment options, prognoses, and recurrence patterns (5).

The composition of immune cells within the BRCA microenvironment plays a crucial role in disease progression and clinical outcome. Innate immune cells (macrophages, neutrophils, dendritic cells, innate lymphoid cells, myeloid-derived suppressor cells, and natural killer cells) and adaptive immune cells (T cells and B cells) within this microenvironment contribute cooperatively to tumor progression and therapeutic response (6). Therefore, elucidating the intricate interplay between molecular characteristics of BRCA and multiple immune cells within the tumor microenvironment can lead to favorable outcomes.

High temperature requirement A (HTRA) proteases are a family of serine proteases involved in various cellular processes (7), among which HTRA1 is the first identified and extensively studied member (8). It plays a crucial role in several fundamental physiological and pathological processes, such as cell proliferation and invasion, and cell adhesion (9). Current studies on HTRA1 generally support its role as a tumor suppressor. Decreased HTRA1 expression is associated with poor prognosis in ovarian cancer (10), metastatic melanoma (11), prostate cancer (12), lung cancer (13), endometrial cancer (14), and BRCA (15). However, the mechanisms underlying the effects of decreased HTRA1 expression and its role in tumor immune-cell infiltration, immune checkpoints, diagnosis, and prognosis of BRCA have not been fully elucidated. In this study, we conducted a comprehensive bioinformatics analysis utilizing The Cancer Genome Atlas (TCGA) (https://portal.gdc.cancer.gov/) and Gene Expression Omnibus (GEO) (https://www.ncbi.nlm.nih.gov/geo/) to assess the important effects of HTRA1 in BRCA (16). We present this article in accordance with the STREGA reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-23-773/rc).

Methods

The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).

Data collection

The RNAseq data in the Fragments Per Kilobase per Million format (level 3 HTSeq-FPKM) and the corresponding clinicopathological information of 1,222 samples (1,109 BRCA and 113 normal breast tissue samples, August 2021) were downloaded from the Breast Invasive Carcinoma Project in TCGA. Then the FPKM format data were converted into transcripts per million (TPM) reads for subsequent analysis. Two BRCA datasets, GSE109169 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE109169) (17) and GSE10780 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE10780) (18) from GEO database were also downloaded for validation of HTRA1 expression. GSE109169 was used to assess the gene expression profiles of 25 paired normal breast and tumor tissues based on the detection platform GPL5175 (HuEx-1_0-st) Affymetrix Human Exon 1.0 ST Array [transcript (gene) version]. GSE10780 was used to analyze the gene expression profiles of 143 normal breast tissues and 42 BRCA tissues based on the detection platform GPL570 (HG-U133_Plus_2) Affymetrix Human Genome U133 Plus 2.0 Array.

Analysis of HTRA1 mRNA expression level in BRCA and normal breast tissue samples

The expression level of HTRA1 mRNA was firstly compared with data from 33 types of human cancer extracted from TCGA (19). For TCGA-BRCA datasets, we further analyzed HTRA1 mRNA expression in tumor samples and their paired normal controls based on BRCA subtypes. The data of HTRA1 extracted from GSE109169 and GSE10780 datasets were used to validate TCGA-BRCA results. Then the TCGA-BRCA samples were divided into high- and low-HTRA1 expression group based on the median expression of HTRA1. Differentially expressed genes (DEGs) between the two groups were screened out and analyzed by “DESeq2” (v1.26.0) R package (20) with log fold change absolute value >1.5 and P value <0.05 as a threshold. The volcano plot and heat map of the DEGs were visualized using the “ggplot2” (v3.3.3) R package.

Functional enrichment analysis of HTRA1-associated DEGs in BRCA

The “ClusterProfiler” (v3.14.3) R package (21) was used for the functional annotation and Gene Set Enrichment Analysis (GSEA) (22) for the selected DEGs. The curated reference gene sets c2. cp.v7.2. symbols.gmt from MSigDB (http://www.gsea-msigdb.org/gsea/msigdb) were designated for GSEA (23). We identified significantly enriched clusters using false discovery rate (FDR) <0.25 and P_adjust<0.05 as a threshold. A protein-protein interaction (PPI) analysis was performed using the STRING database (24) and visualized by Cytoscape software (v3.9.1) (25). The interaction threshold was set to 0.4.

Correlation analysis between HTRA1 expression, immune-cell infiltration, and immuno-related biomarkers in BRCA

The ssGSEA algorithm in the “GSVA” (v1.34.0) R package (26) was used to evaluate the tumor infiltration status of 24 types of immune cell (27). Spearman correlation analysis was performed to determine the relationship between HTRA1 expression, immune-cell infiltration status, and immunity-related biomarkers. Among them, PDCD1 [programmed cell death protein 1 (PD-1)], CD274 [programmed cell death ligand 1 (PD-L1)], CTLA4 and LAG3 are important immune checkpoint markers that are associated with tumor immune escape (28); TIGIT, HAVCR2, LAYN, and CXCL3 are T cell exhaustion markers; IFNG [interferon γ (IFN-γ)], IL2 [interleukin 2 (IL-2)], GZMB (granzyme B), and PRF1 (perforin) are T cell effector function markers. TP53 is a tumor suppressor (29). MKI67 (Ki67) is a proliferation marker.

Analysis of HTRA1 genovariation and DNA methylation in the CpG islands of HTRA1 DNA in BRCA samples

The genomic variation in HTRA1 gene was analyzed by cBioPortal database (https://www.cbioportal.org/) (30,31) including eight BRCA datasets: TCGA (19), METABRIC (32), INSERM (33), MBC project (https://mbcproject.org/), BCCRC (34), Broad (35), MSKCC (36), and Sanger (37). Kaplan-Meier survival analysis and log-rank test were performed to determine the prognostic significance of the genomic variation in HTRA1 gene. The DNA methylation in the CpG sites of HTRA1 gene and their prognostic values were analyzed by MethSurv database (https://biit.cs.ut.ee/methsurv/) (38). The association between CpG methylation of HTRA1 DNA and the overall survival (OS) of BRCA patients was also evaluated.

Correlation analysis between HTRA1 expression and the clinicopathological characteristics of BRCA patients

The clinicopathological information of BRCA patients was extracted from TCGA-BRCA (39). The differences of various clinicopathological parameters were compared between the high- and low-HTRA1 expression group by the Shapiro-Wilk normality test, Kruskal-Wallis test, and Dunn’s multiple hypothesis test. The results were corrected by the Bonferroni method and visualized using the “ggplot2” (v3.3.3) R package. Logistic regression analysis was used to evaluate the relationship between HTRA1 expression and the clinicopathological characteristics of BRCA patients.

Evaluation of the prognostic significance of HTRA1 expression in BRCA patients and their therapeutic response

The survival data of BRCA patients from TCGA-BRCA (39) and GEO (16,40) were analyzed by the “survival” (v3.2-10) and “survminer” (v.0.4.9) R package. Kaplan-Meier survival analysis as well as univariate and multivariate Cox regression were conducted to determine the survival of BRCA patients based on HTRA1 expression. The diagnostic and time-dependent survival ROC curve were generated using the “pROC” (v1.17.0.1) and “timeROC” (v0.4) R package. The online tool ROC Plotter (http://www.rocplot.org/site/index) was used to evaluate the impact of HTRA1 expression on the therapeutic response of breast patients from GEO (41).

Statistical analysis

The differences in HTRA1 expression between groups were analyzed using the Wilcoxon signed-rank test and one-way ANOVA. The correlations between HTRA1 expression and clinicopathological characteristics of patients were assessed using Fisher’s exact test, Chi-squared test, Wilcoxon signed-rank test, and logistic regression. Other statistical analyses were detailed above. P<0.05 was considered to be statistically significant.

Results

HTRA1 expression was significantly dysregulated in multiple types of cancer including BRCA

Figure 1A shows HTRA1 expression in 33 types of cancer in TCGA. The mRNA level of HTRA1 was found to be downregulated in seven types of cancer (BLCA: P=0.003; BRCA: P=0.003, shown also in Figure 1B; CESC: P=0.004; KICH: P<0.001; KIRP: P<0.001; PRAD: P<0.001, and UCEC: P<0.001), and upregulated in eight types of cancer (CHOL: P=0.003; ESCA: P=0.002; HNSC: P<0.001; KIRC: P<0.001; LIHC: P<0.001; PCPG: P=0.015; STAD: P<0.001, and THCA: P<0.001). HTRA1 mRNA expression in breast tumor tissues was also assessed using data from GSE109169 and GSE10780 datasets (Figure 1C,1D, respectively). The comparative results showed that HTRA1 mRNA expression in tumor tissues was significantly lower than that in normal breast tissues (P_adj=0.0233, P_adj=0.0111).

Figure 1 Expression of HTRA1 is significantly dysregulated in BRCA and other cancer types. (A) HTRA1 expression in 33 types of cancer compare with their corresponding adjacent normal tissues. ns, P≥0.05; *, P<0.05; **, P<0.01; ***, P<0.001. Wilcoxon rank sum test. (B) HTRA1 expression was significantly lower in the BRCA tissues compared to the adjacent peritumoral normal breast tissues in TCGA-BRCA. (C,D) HTRA1 expression was significantly lower in breast tumor tissues as compared to the normal breast tissues in (C) GSE109169 (paired samples) and (D) GSE10780. (E-H) Differential expression of HTRA1 in BRCA tissues and matched normal samples from TCGA-BRCA by subtype: Basal, P<0.001; HER2, P=0.4368; Luminal A, P<0.001, and Luminal B, P=0.1448. Paired-sample t-test, two-tailed. (I) Volcano plot of DEGs associated with HTRA1 expression. Every dot corresponds to a gene. The red and dark blue dots represent significant DEGs, and the red dots indicate upregulated expression, dark blue dots indicate downregulated expression, and black dots indicate DEGs with no significant difference. (J) Heat map (top 10 DEGs) showed the level of DEGs in TCGA-BRCA with high- and low-HTRA1 expression divided by the median HTRA1 level. HTRA1, high temperature requirement A1; HER2, human epidermal growth factor receptor 2; ACC, adrenocortical cancer; BLCA, bladder urothelial carcinoma; BRCA, breast cancer; CESC, cervical squamous cell carcinoma and endocervical adenocarcinoma; CHOL, cholangiocarcinoma; COAD, colon adenocarcinoma; DLBC, lymphoid neoplasm diffuse large B-cell lymphoma; ESCA, esophageal carcinoma; GBM, glioblastoma multiforme; HNSC, head and neck squamous cell carcinoma; KICH, kidney chromophobe; KIRC, kidney renal clear cell carcinoma; KIRP, kidney renal papillary cell carcinoma, LAML, acute myeloid leukemia; LGG, lower grade glioma, LIHC, Liver hepatocellular carcinoma; LUAD, lung adenocarcinoma; LUSC, lung squamous cell carcinoma, MESO, mesothelioma; OV, ovarian serous cystadenocarcinoma; PAAD, pancreatic adenocarcinoma; PCPG, pheochromocytoma and paraganglioma; PRAD, Prostate adenocarcinoma; READ, rectal adenocarcinoma; SARC, sarcoma, SKCM, skin cutaneous melanoma; STAD, stomach adenocarcinoma; TGCT, testicular germ cell tumors; THCA, thyroid cancer; THYM, thymoma; UCEC, uterine corpus endometrioid cancer; UCS, uterine carcinosarcoma; UVM, uveal melanoma.

Figure 1E-1H shows HTRA1 expression in tumor tissues and their matched normal controls from TCGA-BRCA by subtype. HTRA1 expression was significantly lower in Basal-like BRCA than that in the paired normal tissues (n=16), while it was not significantly different for either HER2 positive (n=9) or luminal B BRCA (n=21), but was significantly higher for luminal A BRCA (n=56).

The median HTRA1 expression value was used to classify TCGA-BRCA samples into high- and low-HTRA1 expression group. Using the absolute log fold change >1.5 and P<0.05 as a threshold, we screened out 607 DEGs (74 upregulated and 533 downregulated) in the high-HTRA1 expression group compared to the low-HTRA1 expression group (Figure 1I). The ten most significant DEGs are shown in the single gene co-expression heat map in Figure 1J.

Furthermore, the immunohistochemical results from Human Protein Atlas (https://www.proteinatlas.org/) (42,43) verified the expression of HTRA1 protein. In 17 out of 23 samples, HTRA1 staining was absent in breast tumor cells, whereas its weak staining was observed in the cytoplasm or cytomembrane of tumor cells in the remaining six samples (Tables S1,S2).

HTRA1-associated DEGs in BRCA were mainly involved in cell proliferation and cell cycle process

The functional enrichment annotation for the HTRA1-associated DEGs in BRCA was conducted by the “clusterProfiler” R package. The results of Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment, encompassing highly enriched biological process, cellular component, and molecular function (P<0.05), are displayed in Figure 2A and Table S3. The most enriched biological process were extracellular structure organization, extracellular matrix organization, and antimicrobial humoral response. Collagen-containing extracellular matrix, nucleosome, and collagen trimer emerged as the most enriched cellular component. Regarding molecular function, ion channel activity, extracellular matrix structural constituent, and collagen binding exhibited the highest levels of enrichment. The top-ranked KEGG pathways included neuroactive ligand-receptor interaction, alcoholism, and protein digestion and absorption. The PPI network of the DEGs revealed an association between HTRA1 and other genes (Figure 2B and Figure S1). GSEA analysis demonstrated that HTRA1-associated DEGs were significantly enriched in cell proliferation-related clusters (Figure 2C-2I), such as M phase (NES =−1.600, P_adj=0.036, FDR =0.032), G2/M checkpoints (NES =−1.729, P_adj=0.036, FDR =0.032), G0 and early G1 (NES =−1.733, P_adj=0.037, FDR =0.033), cell cycle checkpoints (NES =−1.737, P_adj=0.036, FDR =0.032), cell cycle (NES =−1.625, P_adj=0.036, FDR =0.032; NES =−1.640, P_adj=0.036, FDR =0.032), and DNA replication (NES =−1.731, P_adj=0.036, FDR =0.032). HTRA1-associated DEGs were also enriched in cellular senescence (NES =−1.803, P_adj=0.036, FDR =0.032), DNA double strand break repair (NES =−1.685, P_adj=0.036, FDR =0.032), estrogen dependent gene expression (NES =−1.627, P_adj=0.036, FDR =0.032), DNA double strand break response (NES =−1.804, P_adj=0.036, FDR =0.032), and epigenetic regulation of gene expression (NES =−1.805, P_adj=0.036, FDR =0.032) (Figure 2J-2N).

Figure 2 Functional enrichment results of the HTRA1-associated DEGs in BRCA. (A) GO and KEGG enrichment analysis of HTRA1-associated DEGs showed the enriched BP, CC, and MF. (B) The PPI network of the DEGs. Lines represent protein-protein interaction. The larger the circle, the more the interaction between genes. Red, in the high-HTRA1 expression group, log fold change value >1.5. Green, in the low-HTRA1 expression group, log fold change value <−1.5. (C-N) GSEA analysis showed the altered signaling pathway in BRCA based on HTRA1-associated DEGs. HTRA1, high temperature requirement A1; DEG, differentially expressed gene; BRCA, breast cancer; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; BP, biological process; CC, cellular component; MF, molecular function; PPI, protein-protein interaction; GSEA, gene set enrichment analysis.

HTRA1 expression was correlated with the infiltration of multiple immune cells in BRCA tissues

HTRA1 showed a negative correlation with aDC (R=−0.114, P<0.001) and Th2 cells (R=−0.228, P<0.001), but a positive correlation with CD8 T cells (R=0.201, P<0.001), cytotoxic cells (R=0.091, P=0.003), DC (R=0.137, P<0.001), eosinophils (R=0.267, P<0.001), iDC (R=0.485, P<0.001), macrophages (R=0.316, P<0.001), mast cells (R=0.590, P<0.001), neutrophils (R=0.333, P<0.001), NK CD56bright cells (R=0.140, P<0.001), NK CD56dim cells (R=0.063, P=0.035), NK cells (R=0.537, P<0.001), pDC (R=0.195, P<0.001), T cells (R=0.091, P=0.002), Tcm (R=0.081, P=0.007), Tem (R=0.191, P<0.001), TFH (R=0.099, P<0.001), Tgd (R=0.169, P<0.001), and Th1 cells (R=0.164, P<0.001) (Figure 3A). The tumor infiltration levels of neutrophils, cytotoxic cells, DC, aDC, NK, and Th2 cells in the high- and low-HTRA1 expression group were consistent with the Spearman analysis results (Figure 3B-3G).

Figure 3 The correlation between immune-cell infiltration and HTRA1 expression in BRCA. (A) Spearman analysis result showed the correlation between the infiltration of 24 types of immune cells and HTRA1 expression in BRCA tissues. (B-G) The infiltration level of (B) neutrophils, (C) cytotoxic cells, (D) NK cells, (E) aDC, (F) Th2 cells, and (G) DC in the high- and low-HTRA1 expression group. *, P<0.05; ***, P<0.001. HTRA1, high temperature requirement A1; BRCA, breast cancer; aDC, activated DC; TReg, regulatory T cells; Th17 cells, type 17 Th cells; NK, natural killer cells; Tcm, T central memory; TFH, T follicular helper; DC, dendritic cells; Th1 cells, type 1 Th cells; Tgd, T gamma delta; Tem, T effector memory; pDC, plasmacytoid DC; iDCs, immature DCs; Th2 cells, type 2 Th cells.

Regarding immune cell markers, HTRA1 expression was positively correlated with specific biomarkers for CD8⁺ T cells (CD8A, CD8B), TFH (BCL6), Th1 (T-bet/TBX21, STAT4), Th2 (CCR3, GATA3, STAT5A, STAT6), Th9 (IRF4, PU.1/SPI1, TGFBR2), Th17 (IL-21R, IL-23R, STAT3), Th22 (AHR, CCR10), Treg (CCR8, FOXP3), M1 macrophage (COX2/PTGS2, INOS/NOS2), M2 macrophage (CD206/MRC1, CD115/CSF1R), and TAM (PDCD1LG2, CD80, CD40); whereas it was negatively correlated with biomarkers for B cells (CD38), Th1 (IL2RB2), and Th17 (IL-17A) in BRCA (Table S4).

HTRA1 expression was correlated with immuno-related markers, TP53, and Ki67 in BRCA

In BRCA tissues, HTRA1 expression showed a positive correlation with the levels of CD274 (PD-L1), TIGIT, HAVCR2, LAYN, IL2 (IL-2), PRF1 (perforin), and TP53. Additionally, it exhibited a negative correlation with the levels of LAG3, GZMB (granzyme B), and MKI67 (Ki-67) (Figure 4).

Figure 4 The correlation between HTRA1 expression and immuno-related markers in BRCA (Spearman correlation). (A-D) Correlation of HTRA1 expression and immune checkpoint marker PDCD1 (PD-1), CD274 (PD-L1), CTLA4 and LAG3. (E-H) Correlation of HTRA1 expression and T cell exhaustion marker TIGIT, HAVCR2, LAYN, and CXCL3. (I-L) Correlation of HTRA1 expression and T cell effector function marker IFNG (IFN-γ), IL2 (IL-2), GZMB (granzyme B), and PRF1 (perforin). (M) Correlation of HTRA1 expression and MKI67 (Ki-67). (N) Correlation of HTRA1 expression and TP53. HTRA1, high temperature requirement A1; BRCA, breast cancer; PDCD1 (PD-1), programmed cell death 1; PD-L1, programmed cell death 1 ligand 1; CTLA4, cytotoxic T-lymphocyte associated protein 4; LAG3, lymphocyte activating 3; TIGIT, T cell immunoreceptor with Ig and ITIM domains; HAVCR2, hepatitis A virus cellular receptor 2; LAYN, layilin; CXCL3, C-X-C motif chemokine ligand 3; IFNG (IFN-γ), interferon gamma; IL-2, interleukin 2; GZMB, granzyme B; PRF1 (perforin), perforin 1; MKI67 (Ki67), marker of proliferation Ki-67; TP53, tumor protein p53.

HTRA1 genovariation had no impact on the survival of BRCA patients

Only 1% genovariation was observed in HTRA1 based on data from eight BRCA datasets, including TCGA (n=1,084), METABRIC (n=2,509), INSERM (n=216), MBC project (n=237), BCCRC (n=65), Broad (n=103), MSKCC (n=70), and Sanger (n=100) (Figure 5A). Kaplan-Meier survival curve and log-rank test results showed that there was no significant difference in OS (P=0.726), progression free survival (P=0.456), relapse free status (P=0.543), disease free survival (P=0.430), and disease-specific survival (DSS) (P=0.400) between patients with or without HTRA1 genovariation (Figure 5B-5F). Further analysis revealed the incidence of HTRA1 mutation was 0.42% in MBC project, 0.46% in INSERM, and 0.37% in TCGA BRCA samples, respectively (Figure 5G). The highest level of HTRA1 amplification was observed in the MBC project dataset (2.53%), followed by METABRIC and TCGA (both 0.74%) and INSERM (0.46%). The mRNA expression level of HTRA1 showed a slight increase in tumors with amplified HTRA1 gene compared to those with diploid status (median, 11.7013 vs. 11.4394) (Figure 5H).

Figure 5 HTRA1 genovariation has no impact on HTRA1 expression in breast cancer. (A) OncoPrint visual summary depicted the HTRA1 genovariation. *, not all samples were profiles (altered / profiled =40/3,989). (B-F) Kaplan-Meier survival curve showed (B) overall survival, (C) progression free survival, (D) relapse free status, (E) disease free survival, and (F) disease-specific survival of BRCA patients with altered or unaltered HTRA1 genovariation. (G) HTRA1 genetic alteration in different datasets deposited in cBioPortal. (H) Comparison of HTRA1 expression in BRCA patients of diploid and amplified HTRA1 in METABRIC (n=2,509, Mann-Whitney t-test, P=0.2057). HTRA1, high temperature requirement A1; BRCA, breast cancer; CNA, copy-number alterations.

HTRA1 DNA methylation affected the prognosis of BRCA patients

In the HTRA1 DNA, there were 13 CpG islands showed aberrant DNA methylation, including cg05723130, cg15868400, cg06777002, cg07863439, cg15719652, cg01145353, cg09553346, cg17875153, cg05426956, cg22701672, cg00701951, cg03974204, and cg23791011 (Figure 6). Moreover, the methylation levels of seven CpG islands, namely cg00701951, cg01145353, cg01962937, cg04768425, cg06474225, cg15719652, and cg26466234, were associated with the prognosis of BRCA (Table 1). Increased methylation levels in these seven CpG islands, particularly cg01145353, cg15719652, and cg00701951, were associated with poorer OS in BRCA patients.

Figure 6 Heat map displayed DNA methylation status in HTRA1 DNA in BRCA tissues. HTRA1, high temperature requirement A1; BRCA, breast cancer.

HTRA1 expressions correlated with multiple clinicopathological characteristics in BRCA

Table 2 displays the association between HTRA1 expression and clinicopathological characteristics in BRCA patients. Significant differences were observed between the high- and low-HTRA1 expression group in terms of clinical T stage, histological type, PR status, ER status, PAM50 subtype, and age. Additionally, HTRA1 expression also varied significantly with age, pathologic stage, histological type, T stage, N stage, ER status, PR status, and PAM50 subtype in BRCA (Figure 7). Logistic regression analysis further confirmed a positive correlation of HTRA1 expression with N stage, age, ER status, PR status, histological type, and PAM50 subtype in BRCA (Table 3).

Figure 7 HTRA1 expression correlated with multiple clinicopathological characteristics of BRCA patients. (A-H) The correlation between HTRA1 expression and age (A), pathologic stage (B), histological type (C), T stage (D), ER (E), PR (F), PAM50 (G), and N stage (H) of BRCA patients. *, P<0.05; **, P<0.01, ***, P<0.001. HTRA1, high temperature requirement A1; BRCA, breast cancer; PR, progesterone receptor; ER, estrogen receptor; HER2, human epidermal growth factor receptor 2; PAM50, prediction analysis of microarray 50; LumA, luminal A; LumB; luminal B.

High HTRA1 expression was associated with favorable prognosis of BRCA patients from GEO datasets, but not TCGA

Kaplan-Meier survival curves demonstrated that the relative HTRA1 expression level (high or low) in TCGA-BRCA patients was not correlated with OS (P=0.776), DSS (P=0.551), or progression-free interval (PFI, P=0.446). However, when examining the GEO database (2017, n=5,143, HTRA1 Affymetrix ID is 201185_at) (40,44), high HTRA1 expression was associated with favorable OS and relapse-free survival (RFS) in BRCA patients (Figure 8), indicating a significantly reduced risk for recurrence and death.

Figure 8 The result of KM-plotter analysis showed the impact of HTRA1 expression on breast patients’ survival. (A) Overall survival (n= 1,402); (B) relapse-free survival (n = 3,951). HTRA1, high temperature requirement A1; KM, Kaplan-Meier; HR, hazard ratio.

Furthermore, our ROC curves indicated that HTRA1 expression may not have diagnostic value in distinguishing tumor from adjacent normal tissue [area under the curve (AUC) =0.572]. Additionally, we used the online tool ROC Plotter to assess the impact of HTRA1 expression on therapeutic response in BRCA patients from GEO (41). The selected datasets were all from clinical trials including GSE16446 (NCT00017095, NCT00336791), GSE41998 (NCT00455533), GSE50948 (ISRCTN86043495), GSE66305 (NCT00429299), and GSE16391 (NCT00004205). The results showed that for the luminal B subtype, patients with lower HTRA1 expression had a higher likelihood of achieving pathological complete response to chemotherapy (Figure 9 and Table S5).

Figure 9 The ROC curves showed the differential HTRA1 expression between breast cancer patients with or without pathological complete response to chemotherapy and anti-HER2 therapy. (A) Chemotherapy; (B) anti-HER2 therapy. HTRA1, high temperature requirement A1; ER, estrogen receptor; HER2, human epidermal growth factor receptor 2; TNBC, triple-negative breast cancer.

Discussion

By evaluating HTRA1 expression in human cancers within TCGA database, we observed dysregulation of HTRA1 in 15 out of 33 types of cancer. Specifically, in BRCA, both TCGA and GEO datasets showed a significant reduction in HTRA1 expression in tumor tissues. However, low HTRA1 expression was associated with an unfavorable prognosis in patients from GEO, but not TCGA. These results are consistent with previous studies. Wang et al. reported a notable decrease or absence of HTRA1 expression in in situ or invasive ductal BRCA tissues compared to its prominent expression in normal breast ductal epithelium tissues (45). Lehner et al. measured HTRA1 mRNA expression in 131 BRCA patients (15) and found significantly higher levels of HTRA1 in lower tumor stage, which correlated with longer OS and disease free survival. Franco et al. evaluated 66 sentinel node positive BRCA cases and discovered a negative correlation between HTRA1 expression and metastasis development as well as local relapse (46). However, a retrospective analysis suggested that six epithelial mesenchymal transition genes including LUM, SFRP4, COL6A3, MMP2, CXCL12, and HTRA1 were continuously expressed at higher levels in samples from patients with metastatic BRCA (47). We also observed a significant correlation between HTRA1 expression in BRCA tissues and various factors such as age, clinical N stage, histological type, PR status, ER status, and PAM50 subtype. However, it is important to note that HTRA1 expression exhibits high heterogeneity among different BRCA subtypes. This finding requires further investigation due to limitations in sample size and variations in sampling sites and procedures.

Subsequently, we screened for 607 HTRA1-related DEGs from TCGA-BRCA samples based on median HTRA1 expression. These DEGs primarily participate in extracellular matrix organization processes and exhibit molecular functions related to these processes. GSEA results demonstrated that the HTRA1-related DEGs were predominantly involved in cell proliferation, cell cycle regulation, cellular senescence, and the estrogen-dependent gene expression clusters. Wang et al. also reported decreased HTRA1 expression in BRCA cell lines (45). HTRA1 siRNA-transfected non-tumorigenic MCF10A normal breast epithelium cells exhibited smaller size, increased growth rate, as well as enhanced migration and invasion capabilities. Protease nexin-1 (PN-1) promotes cancer cell metastasis by remodeling the tumor matrix. Tang et al. demonstrated that upregulated PN-1 in BRCA promotes cell invasion, migration, and stemness through the EGF/EGFR/PKCδ/MEK/ERK/EGR1 axis, with HTRA1 being a target of PN-1. Consequently, HTRA1 overexpression inhibits migration and invasion of BRCA cells, while PN-1 overexpression nullifies this inhibitory effect (48). Collectively, these findings confirm the inhibitory impact of HTRA1 on cell proliferation in breast cells. Nevertheless, there are discrepancies among studies, and the current evidence does not definitively establish HTRA1 as a tumor suppressor (9).

Genetic variation has been implicated in the diagnosis, treatment, and prognosis of BRCA patients (49). However, in our analyses, HTRA1 genovariation was not found to be correlated with survival of BRCA patients. We hypothesize that this may be due to the minimal impact of genovariation on HTRA1 expression. DNA methylation is an epigenetic mechanism which plays a crucial role in tumorigenesis. Alterations in HTRA1 DNA methylation can affect its expression in BRCA (15). Notably, there is an inverse correlation between mRNA expression and DNA methylation at the transcription start site of the HTRA1 promoter in BRCA cells (45). In our assessment of the relationship between HTRA1 DNA methylation and prognosis of BRCA patients, we observed that increased methylation at three CpG islands (cg01145353, cg15719652, and cg00701951) was associated with poorer OS compared to decreased methylation.

The immune contexture within the breast tumor microenvironment can influence tumor growth and metastasis. A previous study evaluated tumor-infiltrating lymphocytes (TILs) and their relevance to therapeutic response and prognosis in 3,771 patients with primary BRCA undergoing neoadjuvant chemotherapy (50). The results demonstrated that increased TIL density predicted chemotherapy response across all BRCA subtypes. Furthermore, TIL density was associated with a survival benefit in HER2-positive and TNBC patients (50). Our study revealed potential correlations between HTRA1 expression and the proportion of tumor-infiltrating cells, as well as immune cell markers such as CD8A, CCR3, GATA3, and STAT6. Additionally, significant associations were found between HTRA1 expression and the immune checkpoint PD-L1, T cell exhaustion marker, T cell effector function marker, TP53, and Ki67. These findings support the potential suppressive role of HTRA1 in BRCA. Recent studies have investigated the use of anti-PD-1/PD-L1 agents in BRCA, similar to melanoma or non-small cell lung cancer, particularly in TNBC. Encouraging outcomes have been observed with monotherapy or combination treatment (51). The signature of exhausted T cells can also predict immunotherapy response in BRCA (52). Therefore, our results suggest that HTRA1 may serve as a potential target to enhance the efficacy of immunotherapy in BRCA patients.

Folgueira et al. identified predictive markers for neoadjuvant chemotherapy response in BRCA patients (53). They highlighted a trio of genes, including HTRA1, CLPTM1, and MTSS1 genes, which could accurately differentiate responsive from non-responsive tumors. Another study involving 333 non-metastatic locally advanced BRCA patients who received neoadjuvant chemotherapy demonstrated that higher HTRA1 expression was associated with increased recurrence risk and cancer-related death (54). Our findings indicate that HTRA1 expression can only discriminate between responsive and non-responsive Luminal B patients undergoing chemotherapy. However, the utility of HTRA1 as a biomarker for BRCA prognosis and treatment response requires further comprehensive understanding.

There were some limitations in our study. Firstly, we solely analyzed HTRA1 expression in BRCA samples from TCGA and GEO public database, and it would be beneficial to validate our findings using clinical samples. Secondly, we did not directly assess the molecular mechanism and activity of HTRA1 through in vitro and in vivo experiments, therefore our hypothesis about the effects of genetic variation, expression levels, and progression/prognosis/response to therapy are speculative. Thirdly, the relationships between HTRA1 expression and the function of tumor-infiltrating immune cells within the breast tumor microenvironment needs to be elucidated through well-designed studies. Further investigations are necessary to unravel the complex regulatory role of HTRA1 in BRCA.

Conclusions

In conclusion, our analysis revealed a significant decrease in HTRA1 expression in BRCA datasets. Our findings suggest that HTRA1 may exert its influence on BRCA development through the modulation of genes involved in tumor cell proliferation and immunity. Additionally, we observed that HTRA1 expression might play a crucial role in chemotherapy response and DNA methylation. Thus, HTRA1 expression is associated with BRCA prognosis. Collectively, these findings propose HTRA1 as a potential therapeutic target and a useful prognostic biomarker in BRCA.

Acknowledgments

Funding: This work was supported by the Jilin Province Health Science and Technology Capability Improvement Project (grant No. 2021LC127 to D.Z.).

Footnote

Reporting Checklist: The authors have completed the STREGA reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-23-773/rc

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-23-773/prf

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-23-773/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Sung H, Ferlay J, Siegel RL, et al. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin 2021;71:209-49. [Crossref] [PubMed]
2. Sørlie T, Perou CM, Tibshirani R, et al. Gene expression patterns of breast carcinomas distinguish tumor subclasses with clinical implications. Proc Natl Acad Sci U S A 2001;98:10869-74. [Crossref] [PubMed]
3. Ye F, Dewanjee S, Li Y, et al. Advancements in clinical aspects of targeted therapy and immunotherapy in breast cancer. Mol Cancer 2023;22:105. [Crossref] [PubMed]
4. Huppert LA, Gumusay O, Idossa D, et al. Systemic therapy for hormone receptor-positive/human epidermal growth factor receptor 2-negative early stage and metastatic breast cancer. CA Cancer J Clin 2023;73:480-515. [Crossref] [PubMed]
5. Szymiczek A, Lone A, Akbari MR. Molecular intrinsic versus clinical subtyping in breast cancer: A comprehensive review. Clin Genet 2021;99:613-37. [Crossref] [PubMed]
6. Qayoom H, Sofi S, Mir MA. Targeting tumor microenvironment using tumor-infiltrating lymphocytes as therapeutics against tumorigenesis. Immunol Res 2023;71:588-99. [Crossref] [PubMed]
7. Wang Y, Nie G. Overview of Human HtrA Family Proteases and Their Distinctive Physiological Roles and Unique Involvement in Diseases, Especially Cancer and Pregnancy Complications. Int J Mol Sci 2021;22:10756. [Crossref] [PubMed]
8. Zumbrunn J, Trueb B. Primary structure of a putative serine protease specific for IGF-binding proteins. FEBS Lett 1996;398:187-92. [Crossref] [PubMed]
9. Altobelli E, Marzioni D, Lattanzi A, et al. HtrA1: Its future potential as a novel biomarker for cancer. Oncol Rep 2015;34:555-66. [Crossref] [PubMed]
10. Chien J, Staub J, Hu SI, et al. A candidate tumor suppressor HtrA1 is downregulated in ovarian cancer. Oncogene 2004;23:1636-44. [Crossref] [PubMed]
11. Baldi A, De Luca A, Morini M, et al. The HtrA1 serine protease is down-regulated during human melanoma progression and represses growth of metastatic melanoma cells. Oncogene 2002;21:6684-8. [Crossref] [PubMed]
12. Welsh JB, Sapinoso LM, Su AI, et al. Analysis of gene expression identifies candidate markers and pharmacological targets in prostate cancer. Cancer Res 2001;61:5974-8. [PubMed]
13. Esposito V, Campioni M, De Luca A, et al. Analysis of HtrA1 serine protease expression in human lung cancer. Anticancer Res 2006;26:3455-9. [PubMed]
14. Mullany SA, Moslemi-Kebria M, Rattan R, et al. Expression and functional significance of HtrA1 loss in endometrial cancer. Clin Cancer Res 2011;17:427-36. [Crossref] [PubMed]
15. Lehner A, Magdolen V, Schuster T, et al. Downregulation of serine protease HTRA1 is associated with poor survival in breast cancer. PLoS One 2013;8:e60359. [Crossref] [PubMed]
16. Barrett T, Wilhite SE, Ledoux P, et al. NCBI GEO: archive for functional genomics data sets--update. Nucleic Acids Res 2013;41:D991-5. [Crossref] [PubMed]
17. Chang JW, Kuo WH, Lin CM, et al. Wild-type p53 upregulates an early onset breast cancer-associated gene GAS7 to suppress metastasis via GAS7-CYFIP1-mediated signaling pathway. Oncogene 2018;37:4137-50. [Crossref] [PubMed]
18. Chen DT, Nasir A, Culhane A, et al. Proliferative genes dominate malignancy-risk gene signature in histologically-normal breast tissue. Breast Cancer Res Treat 2010;119:335-46. [Crossref] [PubMed]
19. Hoadley KA, Yau C, Hinoue T, et al. Cell-of-Origin Patterns Dominate the Molecular Classification of 10,000 Tumors from 33 Types of Cancer. Cell 2018;173:291-304.e6. [Crossref] [PubMed]
20. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 2014;15:550. [Crossref] [PubMed]
21. Yu G, Wang LG, Han Y, et al. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS 2012;16:284-7. [Crossref] [PubMed]
22. Subramanian A, Tamayo P, Mootha VK, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A 2005;102:15545-50. [Crossref] [PubMed]
23. Liberzon A, Birger C, Thorvaldsdóttir H, et al. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst 2015;1:417-25. [Crossref] [PubMed]
24. Szklarczyk D, Gable AL, Nastou KC, et al. The STRING database in 2021: customizable protein-protein networks, and functional characterization of user-uploaded gene/measurement sets. Nucleic Acids Res 2021;49:D605-12. [Crossref] [PubMed]
25. Shannon P, Markiel A, Ozier O, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res 2003;13:2498-504. [Crossref] [PubMed]
26. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics 2013;14:7. [Crossref] [PubMed]
27. Bindea G, Mlecnik B, Tosolini M, et al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immunity 2013;39:782-95. [Crossref] [PubMed]
28. Ibrahim R, Saleh K, Chahine C, et al. LAG-3 Inhibitors: Novel Immune Checkpoint Inhibitors Changing the Landscape of Immunotherapy. Biomedicines 2023;11:1878. [Crossref] [PubMed]
29. Wang H, Guo M, Wei H, et al. Targeting p53 pathways: mechanisms, structures, and advances in therapy. Signal Transduct Target Ther 2023;8:92. [Crossref] [PubMed]
30. Cerami E, Gao J, Dogrusoz U, et al. The cBio cancer genomics portal: an open platform for exploring multidimensional cancer genomics data. Cancer Discov 2012;2:401-4. [Crossref] [PubMed]
31. Gao J, Aksoy BA, Dogrusoz U, et al. Integrative analysis of complex cancer genomics and clinical profiles using the cBioPortal. Sci Signal 2013;6:pl1. [Crossref] [PubMed]
32. Pereira B, Chin SF, Rueda OM, et al. The somatic mutation profiles of 2,433 breast cancers refines their genomic and transcriptomic landscapes. Nat Commun 2016;7:11479. [Crossref] [PubMed]
33. Lefebvre C, Bachelot T, Filleron T, et al. Mutational Profile of Metastatic Breast Cancers: A Retrospective Analysis. PLoS Med 2016;13:e1002201. [Crossref] [PubMed]
34. Shah SP, Roth A, Goya R, et al. The clonal and mutational evolution spectrum of primary triple-negative breast cancers. Nature 2012;486:395-9. [Crossref] [PubMed]
35. Banerji S, Cibulskis K, Rangel-Escareno C, et al. Sequence analysis of mutations and translocations across breast cancer subtypes. Nature 2012;486:405-9. [Crossref] [PubMed]
36. Nixon MJ, Formisano L, Mayer IA, et al. PIK3CA and MAP3K1 alterations imply luminal A status and are associated with clinical benefit from pan-PI3K inhibitor buparlisib and letrozole in ER+ metastatic breast cancer. NPJ Breast Cancer 2019;5:31. [Crossref] [PubMed]
37. Stephens PJ, Tarpey PS, Davies H, et al. The landscape of cancer genes and mutational processes in breast cancer. Nature 2012;486:400-4. [Crossref] [PubMed]
38. Modhukur V, Iljasenko T, Metsalu T, et al. MethSurv: a web tool to perform multivariable survival analysis using DNA methylation data. Epigenomics 2018;10:277-88. [Crossref] [PubMed]
39. Liu J, Lichtenberg T, Hoadley KA, et al. An Integrated TCGA Pan-Cancer Clinical Data Resource to Drive High-Quality Survival Outcome Analytics. Cell 2018;173:400-416.e11. [Crossref] [PubMed]
40. Lánczky A, Győrffy B. Web-Based Survival Analysis Tool Tailored for Medical Research (KMplot): Development and Implementation. J Med Internet Res 2021;23:e27633. [Crossref] [PubMed]
41. Fekete JT, Győrffy B. ROCplot.org: Validating predictive biomarkers of chemotherapy/hormonal therapy/anti-HER2 therapy using transcriptomic data of 3,104 breast cancer patients. Int J Cancer 2019;145:3140-51. [Crossref] [PubMed]
42. Uhlén M, Björling E, Agaton C, et al. A human protein atlas for normal and cancer tissues based on antibody proteomics. Mol Cell Proteomics 2005;4:1920-32. [Crossref] [PubMed]
43. Uhlen M, Zhang C, Lee S, et al. A pathology atlas of the human cancer transcriptome. Science 2017;357:eaan2507. [Crossref] [PubMed]
44. Győrffy B. Survival analysis across the entire transcriptome identifies biomarkers with the highest prognostic power in breast cancer. Comput Struct Biotechnol J 2021;19:4101-9. [Crossref] [PubMed]
45. Wang N, Eckert KA, Zomorrodi AR, et al. Down-regulation of HtrA1 activates the epithelial-mesenchymal transition and ATM DNA damage response pathways. PLoS One 2012;7:e39446. [Crossref] [PubMed]
46. Franco R, Collina F, Di Bonito M, et al. HtrA1 loss is related to aggressive behavior parameters in sentinel node positive breast cancer. Histol Histopathol 2015;30:707-14. [PubMed]
47. Wei LY, Zhang XJ, Wang L, et al. A Six-Epithelial-Mesenchymal Transition Gene Signature May Predict Metastasis of Triple-Negative Breast Cancer. Onco Targets Ther 2020;13:6497-509. [Crossref] [PubMed]
48. Tang T, Zhu Q, Li X, et al. Protease Nexin I is a feedback regulator of EGF/PKC/MAPK/EGR1 signaling in breast cancer cells metastasis and stemness. Cell Death Dis 2019;10:649. [Crossref] [PubMed]
49. Harbeck N, Penault-Llorca F, Cortes J, et al. Breast cancer. Nat Rev Dis Primers 2019;5:66. [Crossref] [PubMed]
50. Denkert C, von Minckwitz G, Darb-Esfahani S, et al. Tumour-infiltrating lymphocytes and prognosis in different subtypes of breast cancer: a pooled analysis of 3771 patients treated with neoadjuvant therapy. Lancet Oncol 2018;19:40-50. [Crossref] [PubMed]
51. Emens LA, Cruz C, Eder JP, et al. Long-term Clinical Outcomes and Biomarker Analyses of Atezolizumab Therapy for Patients With Metastatic Triple-Negative Breast Cancer: A Phase 1 Study. JAMA Oncol 2019;5:74-82. [Crossref] [PubMed]
52. Terranova-Barberio M, Pawlowska N, Dhawan M, et al. Exhausted T cell signature predicts immunotherapy response in ER-positive breast cancer. Nat Commun 2020;11:3584. [Crossref] [PubMed]
53. Folgueira MA, Carraro DM, Brentani H, et al. Gene expression profile associated with response to doxorubicin-based therapy in breast cancer. Clin Cancer Res 2005;11:7434-43. [Crossref] [PubMed]
54. Carrara GFA, Evangelista AF, Scapulatempo-Neto C, et al. Analysis of RPL37A, MTSS1, and HTRA1 expression as potential markers for pathologic complete response and survival. Breast Cancer 2021;28:307-20. [Crossref] [PubMed]