3.1 Identification of influenza A target genes
Six influenza A virus target genes (IATGs) were retrieved from the VThunter database and are listed below: Epidermal growth factor receptor (EGFR), Annexin A5 (ANXA5), Calcium voltage-gated channel subunit alpha1 C (CACNA1C), C-type lectin domain family 4 member M (CLEC4M), Cluster of Differentiation 209 (CD209), and UV Radiation Resistance Associated gene (UVRAG). The following analyses were performed based on the six identified IATGs.
3.2 Significant differences in the expression of IATGs between tumor and normal tissues
To identify the clinical impact of IATGs on cancers, the expression and survival analysis of IATGs were performed. The results showed that IATGs were aberrantly expressed in thirteen solid tumors, including BRCA, COAD, LUSC, BLCA, LUAD, KICH, LIHC, PRAD, KIRP, STAD, KIRC, THCA, and HNSC (p < 0.05). Specifically, the expression levels of IATGs, including CD209, CACNA1C, UVRAG and CLEC4M, were significantly downregulated in expression. On the contrary, EGFR and ANXA5 expression levels were remarkedly elevated in multiple cancers (p < 0.05, Fig. 2A). In addition, the risk ratios (HR) of ANXA5 in LUSC, STAD, BLCA, HNSC, MESO, LIHC, LGG; CLEC4M in LUSC, STAD, THCA, THYM, KIRC; CACNA1C in LUSC, STAD, UVM, KIRP, MESO, OV; EGFR in BLCA, PAAD, SARC; CD209 in UVM as well as UVRAG in UCEC were greater than 1, suggesting that ANXA5, CLEC4M, CACNA1C and EGFR’s aberrant expression is a risk factor for the listed cancer types, respectively. However, CLEC4M in ESCA; CACNA1C in LIHC, LGG, ESCA; EGFR in KIRC, CHOL; CD209 in SKCM, and UVRAG in MESO, LGG had risk ratios less than 1, which means these genes were protective factors for the listed cancer types. Overall, the aberrant expression of IATGs may affect the tumorigenesis and prognosis of diverse cancer types (Fig. 2B).
3.3 Somatic mutations in IATGs
All cancer cells harbor somatic mutations. A certain number of these somatic alterations, generally known as “driver” mutations, confer selective cell subgroup growth advantages and are causally linked to oncogenesis [30, 31]. Thus, the SNV data associated with IATGs was further analyzed to observe the mutation frequency and type in each cancer. The data showed that somatic mutations of IATGs existed in most of the cancer types except UCS, UVM, KICH, PCPG, THCA, CHOL, LAMI, and MESO. Among them, CACNA1C in SKCM as well as EGFR in GBM had the highest SNV frequency up to 20%, which was presented with the darkest labeling color (Fig. 3A). Following that, the missense mutations which were the predominant type of somatic mutations was focused. In detail, SNV percentage analysis revealed the following gene mutation rates: CACNA1C (46%), EGFR (40%), UVRAG (11%), CD209 (9%), CLEC4M (9%) and ANXA5 (7%). The high frequency of mutations in IATGs was observed in these cancer types, including STAD, SKCM, GBM, PAAD, and LUSC (Fig. 3B). In addition, the HRs of UVRAG, CACNA1C, and EGFR in UCE were less than 1, which indicated that they were protective factors for the cancer types. However, the HRs of EGFR and CD209 in LGG, EGFR in COAD, and CACNA1C in KIRP were greater than 1, which suggested that they are risk factors for certain cancer types, respectively. Collectively, these findings revealed a substantial difference in survival between mutant and non-mutated IATG genes, which may impact the prognosis of cancer patients (Fig. 3C).
3.4 Copy number variation of IATGs
Copy number variation (CNV) refers to a phenomenon in which the number of copies of a specific DNA segment varies between individuals. Changes in copy number are tightly associated with cancer initiation and progression [32]. In order to determine the changes in CNV, the CNV data was therefore analyzed for IATGs. It showed that the major types of CNV in IATGs were heterozygous amplification and deletion. The highest percentage of heterozygous amplification is EGFR in TGCT (Fig. 4A).
Additionally, CACNA1C in KIRC, LUSC, and BRCA, and CD209 expression in THCT and ACC were negatively correlated with CNV. However, in most tumors, UVRAG, EGFR, and ANXA5 expression levels were positively correlated with CNV (Fig. 4B). Specifically, EGFR, ANXA5, CACNA1C, CLEC4M, CD209 and UVRAG in UCEC; UVRAG, ANXA5, CACNA1C, CLEC4M, CD209 in KIRC; UVRAG, ANXA5, CLEC4M, CD209 in ACC; EGFR, UVRAG, ANXA5, CLEC4M, CD209 in HNSC; UVRAG, ANXA5 in KIRP; UVRAG, CACNA1C in LAML; UVRAG, ANXA5, CLEC4M, CD209 in MESO; EGFR, UVRAG, ANXA5 in KIRP; ANXA5, CACNA1C in LAML; UVRAG, ANXA5 in MESO; EGFR, UVRAG in CHOL, LGG, LUAD and LUSC; EGFR in GBM; CACNA1C in KICH; UVRAG in PAAD; CACNA1C in PCPG all had a HR which was greater than 1 and defined as risk factors. In contrast, EGFR in LAML and MESO; CLEC4M, CD209 in CESC; ANXA5 in COAD and THCA; and CACNA1C in THYM had a HR less than 1 and were defined as protective factors. Together, these findings illustrate that the vast majority of CNV changes in IATGs are heterozygous amplifications as well as deletions, which correlated with the tumor progress and prognosis (Fig. 4C).
3.5 Methylation analysis of IATGs
DNA methylation is characterized as a key regulator of gene expression. Aberrant DNA methylation profiles are generally linked with cancer events [33, 34]. Therefore, we further investigated the methylation status of IATGs in various cancer types. Firstly, the mRNA expressions of most IATGs were negatively correlated with their methylation levels in most tumors, while CD209 in THYM, HNSC and SKCM, as well as CLEC4M in SKCM, were positively correlated with their methylation levels (p < 0.05, Fig. 5A). Secondly, the hypermethylation of UVRAG in MESO and EGFR in KIRC was associated with low survival, which had risk factors greater than 1 (HR > 1). On the contrary, hypermethylation of UVRAG in LGG; EGFR in LGG and BLCA; CLEC4M in SARC, STAD, and MESO; CD209 in LGG, UVM, SARC, and SKCM; CACNA1C in LGG, UVM, KIRP, and STAD; ANXA5 in LGG, BLCA, LIHC, SKCM, LUAD, GBM, KIRC, and HNSC were associated with high survival and defined as protective factors (HR < 1, Fig. 5B). Lastly, the methylation of IATGs in different cancer types was highly heterogeneous. More hypermethylated genes than hypomethylated genes were found in BRCA, UCEC and LUAD. In contrast, HNSC, ESCA, PAAD, THCA, LUSC, LIHC, and KIRC cancer types had more hypomethylated genes (Fig. 6A).
3.6 Analysis of the relationship between IATGs and cancer-related pathways
Here, the relationship between IATGs and cancer-related pathways were analyzed. Our results revealed that five IATGs included in this study except CLEC4M were significantly involved in the classic cancer-related signaling pathways. The detailed signaling pathways are listed as follows: apoptosis, cell cycle, DNA damage, EMT, hormone AR, hormone ER, PI3K/AKT, RAS/MAPK, RTK, and TSC/mTOR. The numbers in each cell in panel B indicated the percentage of a certain cancer type where the gene was associated with a specific pathway. The major affected pathways for UVRAG were DNA Damage (19%) and EMT (16%) pathway activation. The major involved pathways for EGFR were RTK activation (31%) and DNA Damage inhibition (28%); for CD209, the affected pathways were apoptosis (28%) and EMT (34%) activation; for CACNA1C, the main pathways were apoptosis (31%) and cell cycle (41%) inhibition; for ANXA5, the main pathways involved were EMT activation (19%) and Hormone AR inhibition (19%). These results suggest that IATGs play vital roles in the regulation of cancer-related signaling pathways (Fig. 6B).
3.7 IATGs and tumor drug resistance analysis
Cancer drug resistance remains the biggest challenge in cancer therapy. Overcoming the resistance is still extremely urgent. Here, two chemotherapy drugs (Tanespimycin and Docetaxe), which can be applied to several cancer types, were negatively correlated with the expression of EGFR and ANXA5, but positively correlated with the expression of CACNA1C, UVRAG, CLEC4M, and CD209. Furthermore, the resistance of the following drugs, including BHG712, THZ-2-102-1, TL-1-85, TPCA-1, ZSTK474, AR-42, AT-7519, BMS345541, BX-912, CAY10603, CP466722, I-BET-762, JW-7-24-1, KIN001-102, KIN001-260, Methotrexate, NG-25, NPK76 -II-72-1, Navitoclax, OSI-027, PHA-793887, PI-103, PIK-93, QL-XI-92, TG101348, Tubastatin A, Vorinostat, and WZ3105 were negatively correlated with EGFR, ANXA5 expression and positively correlated with UVRAG, CLEC4M and CD209 expression, respectively. These results suggest that aberrant expression of IATGs may deeply participate in tumor resistance to chemotherapy as well as targeted drug therapy (Fig. 6C).
3.8 Summary of gene set enrichment scores in IATGs tumors
The enrichment score (ES) obtained from the gene set enrichment analysis (GSEA) reflects the degree of specific gene enrichment. The ES of IATGs was upregulated in HNSC, THCA, and ESCA, which were more likely to be enriched in these three cancer types, but downregulated in LUSC, BLCA, PRAD, KIRP, LUAD, COAD, LIHC, and BRCA, which were considered to have a lower possibility of enrichment in these cancer types (Fig. 6D).
3.9 Relationship between IATGs and immunization
At the current stage, cancer immunotherapy has been widely applied in the treatment of diverse cancer types with a good prognosis [35, 36]. Thus, exploring the links between the IATGs and immunity is meaningful. As shown in Fig. 7A, based on the GSVA scores of IATGs, it was indicated that IATGs were positively correlated with the infiltration score in most of the tumor types as well as with the expression of markers of immune infiltrating cells, such as CD4+ T cells, natural killer cells (NK), natural killer T cells (NKT), T follicular helper cells (Tfh), central memory T cells (TCM cells), dendritic cells (DCs) and macrophages. On the contrary, IATGs were negatively correlated with the expression of markers of effector memory T cells (TEM cells), neutrophils, B cell and naive CD4+ T cells (*p value ≤ 0.05; #: FDR ≤ 0.05). Furthermore, correlation analysis of IATGs and the immunostimulatory pathways using the GEPIA2 database revealed a positive correlation between IATGs and the immunostimulatory pathway (Fig. 7B, R = 0.48), the MHC immune pathway (Fig. 7C, R = 0.36), and the chemokine immune pathway (Fig. 7D, R = 0.40). Taken together, these results suggest that IATGs are closely linked with tumor immunity and affect tumor immune pathways to influence tumor progression.
3.10 Confirmation the transcript expression of IATGs in KIRC
In order to increase the reliability of the analyses we performed, we further confirmed the gene expression of ITAGs including EGFR, ANXA5, CACNA1C, CLEC4M, CD209, and UVRAG by qRT-PCR in kidney renal clear cell carcinoma (KIRC). KIRC is one of the most common cancer types and has been chosen as a representative cancer type in this study. HK-2 cell line is a normal kidney epithelial cell line which is used as a control here. 786-O, ACHN and Caki-1 cell lines are the kidney cell carcinoma cell lines. RT-PCR results in the KIRC cell lines showed that EGFR, ANXA5, and CACNA1C were upregulated and CLEC4M and CD209 were downregulated, which were consistent with the data analyzed from the TCGA and GTEx databases (Fig. 8).