An Anoikis-related LncRNA Signature: Integrated Analysis of Immune Infiltration Landscape in Patients with Liver Hepatocellular Carcinoma and Immunotherapy

doi:10.21203/rs.3.rs-3141131/v1

Download PDF

Research Article

An Anoikis-related LncRNA Signature: Integrated Analysis of Immune Infiltration Landscape in Patients with Liver Hepatocellular Carcinoma and Immunotherapy

https://doi.org/10.21203/rs.3.rs-3141131/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Liver hepatocellular carcinoma (LIHC) is one of the most common malignancies worldwide and any factor causing metastasis contributes to its poor prognosis. Long non-coding RNA (lncRNA) plays an important role in promoting the poor progression of LIHC. And the occurrence of anoikis resistance is also essential for tumor metastasis. The study aimed to establish an anoikis-related lncRNA (ARlncRNA) signature for LIHC. The raw data were obtained from the TCGA and GeneCards databases. Using the univariate Cox methods (p < 0.05) and the LASSO (10-fold cross-validation) regression analysis, we first identified a lncRNA-associated signature which consisted of AC100812.1, AL365295.1, AC073352.1, ELFN1-AS1, LINC00513, and MIR4435-2HG. Furthermore, we constructed a risk model based on the ARlncRNAs scores and evaluated it from different perspectives, including survival prognosis, clinical characteristics, signaling pathways, levels of immune cell infiltration, and drug sensitivity. The results of validation indicated that the prognostic performance of the lncRNA-associated risk score is excellent and is more suitable for constructing a prognostic model than other clinical features. We finally identified 3 LIHC subtypes and made predictions on their immune infiltration landscape and drug sensitivity. Our study elucidated the mechanisms of LIHC metastasis and its immune infiltration landscape, which was of great significance in expanding our understanding of LIHC progression and facilitating patients’ personalized management.

liver hepatocellular carcinoma

anoikis

lncRNA

risk model

immune infiltration landscape

immunotherapy

Primary liver cancer is the fifth most common malignancy in the world and ranks third among cancer-related causes of death(Vogel et al. 2022). LIHC represents the most frequent (90%) histological type of primary liver cancer(Lafaro et al. 2015). A data by Dr Christina Fitzmaurice showed that in 2015 alone there were 854,000 newly diagnosed cases of LIHC and 810,000 death cases of LIHC worldwide, meaning that the mortality rate of LIHC at that time was almost equal to its incidence(Collaboration et al. 2017). The situation has now become even worse. The incidence rate of LIHC is still rising and by 2025, the World Health Organization estimates that more than one million people will be diagnosed with LIHC each year(Llovet et al. 2021). Under current conditions, the average 5-year survival rate for LIHC is approximately 18%(Vogel et al. 2022). As a result, LIHC will pose a significant threat to global health in the future. The preferred treatment for LIHC patients with early-stage is surgical resection, followed by liver transplantation, percutaneous ablation and trans-arterial embolization chemotherapy, while advanced patients benefit mainly from immune checkpoint inhibitor-based therapies(Giuliani and Colucci 2009). There is one point deserves additional attention that the prognosis for LIHC patients with early diagnosis and treatment (5-year survival rate: 89%-93%) is significantly better than that of advanced patients with LIHC (survival time < 6 months) (Llovet et al. 2004). The five LIHC stages defined by The Barcelona Clinic Liver Cancer (BCLC) system reflected similar results, with the median survival time for stage-0 was 38 months but less than 6 months for stage-D(Han and Kim 2015; Guiu et al. 2022). Unfortunately, due to the latent onset of LIHC and mild symptoms in the early stage, patients are not diagnosed until the condition is severe enough, resulting in the best opportunity for treatment being lost. The development of novel biomarkers that reflect the severity of the disease is of great significance for patients with LIHC, both for early screening and for later personalized management. The recent boom in gene sequencing technology has created the conditions for this initiative to be realized.

Anoikis (Greek for homelessness) is a programmed pattern of cell death. When normal epithelial cells lose contact with the extracellular matrix (ECM) or adhere to other ECM, anoikis is triggered, through intracellular (perturbation of mitochondria) and extracellular (activation of cell surface death receptors) pathways(Paoli et al. 2013). Anoikis aims to maintain the correct cell position in the corresponding tissues, but for malignant tumors, its triggering failure often occurs, which is called anoikis resistance(Gilmore 2005). Tumors have evolved anoikis resistance through various mechanisms(Taddei et al. 2012). There is evidence that the survival prognosis of tumor cells in the blood, lymphatic circulation and distant organs is directly correlated with anoikis resistance(Zhong and Rescorla 2012). Jia Song et al. confirmed that by suppressing the overexpression of anoikis resistance-associated protein 14-3-3σ, the clinical symptoms and prognosis of LIHC patients could be improved(Song et al. 2021). Similarly found by Suhong Xia et al. that histidine-rich calcium-binding protein (HRC), an anoikis resistance-associated gene, was also overexpressed in LIHC patients and suggested that this promoted the metastasis of LIHC cells in vivo(Xia et al. 2021). Therefore, deciphering the regulatory mechanisms of anoikis resistance would be helpful in improving the clinical management of LIHC.

Long non-coding RNA (lncRNA) is defined as a class of DNA transcription products that are > 200 nt in length and do not have a protein-coding function(Hangauer et al. 2013). Accumulating evidence suggests that by regulating the expression of cancer-related genes at the transcriptional, translational and post-translational levels, lncRNAs can influence a variety of biological activities in malignant tumors, including tumorigenesis, tumor cell proliferation, invasion, EMT, angiogenesis and drug resistance(Peng et al. 2017; Liu et al. 2018; Xu et al. 2018). Given that lncRNAs are not easily interfered with by RNases and their content are stable, we can take non-invasive means to detect the expression levels of extracellular lncRNAs of cancer patients(Shi et al. 2016). Therefore, lncRNAs have great advantages in being novel biomarkers of tumors. To date, a large quantity of lncRNA-related models have been established for many tumors, including thyroid, breast and lung squamous cancers, which is of great significance for revealing the immune infiltrative landscape and providing personalized therapy for patients(Li et al. 2022b; Zhang et al. 2022b; Chen et al. 2022). However, no studies have yet focused on anoikis-related lncRNA models in LIHC.

In the study, we first identified an anoikis-related signature consisting of 6 lncRNAs. Then, we constructed a risk model that was helpful in the prognosis prediction of LIHC. According to the median risk score in the model, we divided all TCGA-LIHC samples into different risk groups. Next, we validated the model and evaluated it from different aspects, including survival prognosis, clinical characteristics, signaling pathways, levels of immune cell infiltration and drug sensitivity. Finally, we identified 3 LIHC subtypes and made predictions on their immune infiltration landscape and drug sensitivity. Through the comprehensive analysis of the genome-wide data in the TCGA database, we constructed an ARlncRNA signature, and further shed light on the different immune infiltration landscapes and personalized immunotherapy schemes of LIHC patients.

Collection and pre-processing of LIHC sample data

The Genomic Data Commons (GDC, https://portal.gdc.cancer.gov) database is a set of cancer data-sharing platforms established by the National Cancer Institute that integrates cancer genomics data from multiple databases, including TCGA. We first captured transcriptomic RNA-seq data from the GDC database for the TCGA-LIHC cohort and its clinical data, and considered raw expression data from para-cancerous tissues as controls. After removing samples with missing or duplicate information, we distinguished samples from the TCGA-LIHC cohort into the control group (n = 50) or the tumor group (n = 374) with the following sample inclusion criteria: 1. patients diagnosed with LIHC; 2. survival time > 30 days; and 3. complete prognostic information. We then performed uniform log2 normalization and ID transformation of the transcriptomic data, and further identified mRNA and lncRNA data in the TCGA-LIHC expression matrix.

Screening of anoikis-related genes and lncRNA data

Anoikis-related genes (ARGs, Relevance scored ≧ 1), obtained from the GeneCards database, were used to contrast with the TCGA-LIHC expression matrix and to output their expression matrix. It is well known that lncRNAs are a class of RNAs that play an important role in protein translation and are intrinsically linked to mRNAs. With a correlation coefficient > 0.4 and p < 0.001, we filtered the ARGs and lncRNA expression matrices by Pearson correlation analysis to synthesize a data matrix of anoikis-related lncRNAs (ARlncRNAs). Using the igraph package, we mapped the regulatory network between ARGs and ARlncRNAs. In addition, we screened for ARlncRNAs that were differentially expressed in the TCGA-LIHC cohort (∣log₂FC∣>1, FDR < 0.05).

Construction of an LIHC-related risk model based on ARlncRNAs scores

We read survival the survival and clinical data of tumor samples from the TCGA-LIHC cohort after removing the normal samples. Using the batch survival univariate Cox methods (p < 0.05), we first screened for ARlncRNAs associated with prognosis in patients with LIHC. Then, we carried out the Logistic Regression with the Least Absolute Shrinkage and Selection Operator (LASSO, 10-fold cross-validation) regression analysis, and further identified an ARlncRNA-related signature for the establishment of the risk model. The formula of risk score (RS) is as follows.

RS=\({\sum }_{k=1}^{n}expr\left({ASlncRNA}^{k}\right)*coef\left({ASlncRNA}^{k}\right)\)

We applied the R language to randomly divide the 374 tumor samples of the TCGA-LIHC cohort into the train and test sets on a 1:1 basis. And all samples were divided into high-risk groups or low-risk groups according to the median risk score of the model. Next, by evaluating the risk score, survival status, and survival prognosis of the 3 cohorts respectively, we verified the accuracy of the risk model in the prognosis prediction of LIHC. The ROC curves plotted by timeROC package were used to evaluate the predicted value of the risk model, and the area under the curve (AUC) is the reference. Based on different sample cohorts, we plotted ROC curves for the risk models at 1-, 3-, and 5-years respectively. Ultimately, we compared the predictive values between risk factors and other clinical characteristics with the help of the ROC curve.

Verification of risk model based on internal data sets

We matched the clinical data of 374 TCGA-LIHC tumor samples with their risk scores, and removed the samples with missing information. The univariate Cox regression analysis identified all factors associated with survival prognosis, and multivariate Cox regression analysis confirmed whether the risk score could be independent of other factors in the prediction of survival prognosis (p < 0.05). Based on the different clinical characteristics of the TCGA-LIHC tumor samples, we further explored the prognostic value of the risk model. The nomogram was plotted to foresee the 1-, 3- and 5-year survival time of LIHC patients by scoring the clinical characteristics and risk factors of the tumor samples. And the calibration curve was used to demonstrate the predictive performance of the nomogram.

GSEA enrichment analysis

The GSEA enrichment analysis can be used to determine whether the different gene sets differ significantly across biological states. Given the risk scores mentioned above, the TCGA-LIHC tumor samples were differentiated into high- and low-risk groups. Using GSEA analysis software, we scored the KEGG signaling pathways differentially expressed between the different risk groups and selected the highest-scoring entries (p < 0.05).

Immune infiltration landscape analysis of the different risk groups

By reading and comparing the information of the risk score file and immune cell infiltration file in TCGA-LIHC tumor samples, we captured the relationships between different risk groups and different immune cell subpopulations. Immune cell subpopulations included XCELL, TIMER, QUANTISEQ, MCPcounter, EPIC, CIBERSORT-ABS, and CIBERSORT. ARlncRNA-associated immunocyte infiltration levels in the TCGA-LIHC cohort were calculated by the GSVA package. Through reading the tumor microenvironment (TME) scoring files of the tumor samples, we analyzed the differences in the abundance of stromal and immune cells between different risk groups, and considered P < 0.05 to be statistically significant. Immune checkpoints, as a class of immunosuppressive molecules, can modulate the immune response and thus influence tumor progression. With the help of ggpubr and ggplot2 packages, we compared the activation of 47 immune checkpoints between the two risk groups.

Drug sensitivity analysis among different risk groups

Half maximal inhibitory concentration (IC50) represents the drug concentration corresponding to 50% death of tumor cells and can be used to predict the degree of drug resistance of tumor cells. Using the pRRophetic and limma packages, we calculated the sensitivity of tumor samples from different risk groups in the TCGA-LIHC cohort to 138 drugs, aiming to identify drugs specific to different risk groups.

Identification of LIHC subtypes

Based on the ConsensusClusterPlus package and ARlncRNA expression matrix, we conducted consensus clustering analysis on the TCGA-LIHC cohort to obtain the best LIHC subtypes that the lncRNAs have the strongest association within the group and the weakest association among the groups. Comparison of survival prognosis between different LIHC isoforms was achieved using the survival and survminer packages. Using the Rtsne and ggplot2 packages, we plotted PCA and T-SNE plots to distinguish the distribution of TGCA-LIHC samples across different LIHC subtypes and different risk groups. The Sankey diagram can be used to demonstrate the potential relationships of TGCA-LIHC samples across different LIHC subtypes and different risk subgroups.

Immune infiltration landscape and drug sensitivity analysis of LIHC subtypes

We analyzed the expression abundance of stromal cells, immune cells, and TME in different LIHC subtypes with reference to the TME scoring file, and observed whether the results were statistically different (p < 0.05). By comparing the information in the LIHC typing results file with the information in the immune cell infiltration file, we obtained the association between different LIHC subtypes and seven immune cell subsets, and further showed the results by a heatmap. The difference in activation of 47 immune checkpoints among different LIHC subtypes was shown by a box plot. We performed drug sensitivity analysis again to predict the drug resistance of patients with different LIHC subtypes to 138 drugs, which is also helpful for the personalized management of LIHC patients.

Statistical analysis

The RNA-seq transcriptome data were preprocessed using the Perl programming language (version 5.32.1) for screening prognosis-related genes. All data visualizations and statistical analyses were performed using R software (version 4.2.0) and a variety of R packages. p < 0.05 was considered statistically significant.

Identification of ARlncRNAs in TCGA-LIHC samples

The TCGA-LIHC gene expression matrix containing 374 tumor samples and 50 para-cancer samples obtained from the GDC database was distinguished into two biological types, mRNA and lncRNA. By comparing 363 ARGs downloaded from the GeneCards database with genes in the mRNA expression matrix, we got 346 ARGs expressed in TCGA-LIHC samples. And through co-expression analysis, we obtained 2902 ARlncRNAs with a regulatory relationship with 346 ARGs (correlation coefficient > 0.4 and p < 0.001). Figure 1A displayed the regulatory relationships between 331 ARGs and 2902 lncRNAs. Among them, 2143 ARlncRNAs were up-regulated and 20 were down-regulated (∣log2FC∣>1, FDR < 0.05) in the TCGA-LIHC cohort (Fig. 1B). The heatmap suggested that differentially expressed ARlncRNAs could well distinguish normal and tumor samples in the TCGA-LIHC cohort (Fig. 1C).

Figure.1 Identification of ARlncRNAs in LIHC. (A) ARGs-ARlncRNAs regulatory network. (B) Volcano plot of ARlncRNAs (∣log₂FC∣>1, FDR < 0.05). (C) Heatmap of the expression of the 50 most significantly up- and down-regulated ARlncRNAs in LIHC.

Identification of prognosis-related ARlncRNAs

After removing normal samples from the TCGA-LIHC cohort, we screened for tumor samples with complete survival data and their differentially expressed lncRNA matrix files. The batch survival univariate Cox method (p < 0.05) was performed on 2163 ARlncRNAs, and we screened 27 ARlncRNAs that had an impact on prognosis (p < 0.001) (Fig. 2A). Sankey diagram demonstrated the regulatory relationships between ARGs and the 27 prognosis-related ARlncRNAs (Fig. 2E). The results showed negative regulatory relationships between ARGs-ARlncRNAs in only 14 pairs (MAPK3 to AC115619.1; RHOA to AC115619.1; HRAS to TMEM220-AS1; AKT1 to RAB11B-AS1; PAK2 to RAB11B-AS1; EEF2K to RAB11B-AS1; PTPN11 to ZSCAN16-AS1; ANXA5 to AC026803.3; ATF4 to AC026803.3; YWHAZ to AC026803.3; PLG to MAPKAPK5-AS1; PLG to CYTOR; PLG to AC011468.1; PLG to AC024060.2). The LASSO regression analysis further identified 6 ARlncRNAs that were used to construct the prognosis model of LIHC, including AC100812.1, AL365295.1, AC073352.1, ELFN1-AS1, LINC00513, and MIR4435-2HG (Fig. 2B-C). The details of ARlncRNAs were demonstrated in Table 1. The heat map suggested that the ARlncRNA signature was significantly down-regulated in the control group, while up-regulated in the tumor group (Fig. 2D).

Table 1

Details of the ARlncRNA signature.
LncRNA	Coefficient	HR	95%HR	p Value
AC100812.1	0.749	3.949	1.896–8.226	< 0.001
AL365295.1	1.454	6.348	2.115–19.051	< 0.001
AC073352.1	0.585	2.916	1.770–4.803	2.65E-05
ELFN1-AS1	0.174	1.356	1.145–1.606	< 0.001
LINC00513	0.573	1.987	1.427–2.765	4.69E-05
MIR4435-2HG	0.564	2.978	1.701–5.216	< 0.001

Construction and internal validation of the LIHC-related risk model associated with ARlncRNAs

After determining the ARlncRNA signature, we calculated risk scores for all samples in 3 cohorts and distinguished them into high- or low-risk groups with reference to the median value of the model risk score (RS).

RS = AC100812.1* 0.749256106692893 + AL365295.1* 1.45427488505416 + AC073352.1* 0.585423143197508 + ELFN1-AS1*0.174483348394598 + LINC00513* 0.572529045824369 + MIR4435-2HG*0.564436019020331.

The results indicated that in all the cohorts, the mortality rate of LIHC patients was higher in the high-risk group (Fig. 3A-B). And the expression of ARlncRNA was significantly increased in the high-risk group (Fig. 3C). It is hinted that the up-regulation of ARlncRNA expression is directly proportional to the poor prognosis of LIHC. K-M plots showed that overall survival (OS) was always significantly lower for patients in the high-risk group than for those in the low-risk group in train, test, and entire sets, implying that factors associated with the high-risk score also contributed to the poor prognosis of LIHC patients (Fig. 3D). The area under the curves (AUCs) can be used to assess the predictive value of the model. In the train set, the AUCs were 0.788 (1-year), 0.808 (3-years), and 0.810 (5-years); the AUCs were 0.684 (1-year), 0.541 (3-years), and 0.596 (5 years) in the test set; as for the entire set, the AUCs were 0.732 (1-year), 0.675 (3-years), and 0.700 (5-years), respectively (Fig. 3E-G). The results indicated that the risk model constructed based on the ARlncRNA signature had a reliable prediction performance in terms of the prognosis of all sets. Figure 3H confirmed that the risk score has a greater prognostic value for patients with LIHC compared to the other clinical features (AUC = 0.732).

Prognostic analysis of risk model in different clinical subgroups

We explored the prognostic value of the risk model across different subtypes of clinical characteristics (gender, age, T, N, M, and stage). The results showed that the survival prognosis of patients in the low-risk group was always better than that of the high-risk group in all clinical characteristics subgroups (Fig. 4). Notably, the prognostic value of our constructed risk model was greater for patients with early-stage LIHC (p < 0.001).

Assessment of the risk model

The results of univariate Cox regression analysis suggested that the prognosis of TCGA-LIHC patients was associated with the stage (p < 0.001; HR = 1.680; 95% CI = 1.369–2.062) and riskScore (p < 0.001; HR = 1.069; 95%CI = 1.053–1.132) (Fig. 5A). The results of multivariate Cox regression analysis provided evidence that riskScore (p = 0.003; HR = 1.064; 95% CI = 1.021–1.109) was an independent prognostic factor for LIHC (Fig. 5B). The curves of the c-index indicated that the riskScore-based model is more accurate in assessing patient prognosis than other clinical characteristics (Fig. 5E). We constructed a nomogram to predict the survival time of LIHC patients at 1-, 3-, and 5-years by scoring age, gender, T, N, M, stage, grade, and risk of the total TCGA-LIHC sample. The results showed that the survival time of LIHC patients at 1-, 3-, and 5-years could reach 0.899, 0.801, and 0.728 respectively (Fig. 5C). The calibration and DCA curves confirmed the value of the nomogram in predicting patient survival time and suggested that its predictive performance is superior to other clinical traits (Fig. 5D and Fig. 5F). In addition, we also confirmed that our identified anoikis-related lncRNA signature had good prognostic performance in predicting the prognosis of LIHC patients compared with the anoikis-related gene signature identified by Yutong Chen et al. and the pyroptosis-related gene signature identified by Min Deng et al. The c-index of our identified anoikis-associated lncRNA signature in predicting the prognosis of LIHC patients was 0.671 (95%CI: 0.559–0.767). The model of Yutong Chen et al., based on the signature of the anoikis-related gene, has a c-index of 0.687 (95% CI: 0.574–0.780); and the c-index of the model constructed by Min Deng et al. based on the signature of pyroptosis-related genes is 0.658(95%CI: 0.545–0.756) (Deng et al. 2022; Chen et al. 2023).

GSEA enrichment analysis

Using GSEA software, we first screened all the Kyoto Encyclopedia of Genes and Genomes (KEGG) signaling pathways involved in ARlncRNAs and ranked their relevance. The results suggested that the KEGG signaling pathways associated with the high-risk group mainly include OOCYTE_MEIOSIS, PYRIMIDINE_METABOLISM, RNA_DEGRADATION, SPLICEOSOME, and CELL_CYCLE (Fig. 6A), while the KEGG signaling pathways involved in the low-risk group mainly include PRIMARY_BILE_ACID_BIOSYNTHESIS, COMPLEMENT_AND_COAGULATION_CASCADES, FATTY_ACID_METABOLISM, DRUG_METABOLISM_CYTOCHROME_P450, and RETINOL _METABOLISM (Fig. 6B).

Immune infiltration landscape analysis and clinical treatment exploration of different ARlncRNAs scored groups

For TCGA-LIHC tumor samples we performed immune cell correlation analysis based on different platform sources (Fig. 7A). The results showed that XCELL platform-derived T cell CD4 + Th2; TIMER platform-derived T cell CD4+, Neutrophil, Macrophage, Myeloid dendritic cells; QUANTISEQ platform-derived Monocyte; EPIC platform-derived uncharacterized cells; CIBERSORT-ABS platform-derived T cell regulatory (Tregs), Macrophage M0 and M2; CIBERSORT platform-derived T cell regulatory (Tregs) and Macrophage M0 all showed a positive regulatory relationship with risk scores. As for the immune cells that showed a negative regulatory relationship with a risk score, the main ones were Endothelial cell, Hematopoietic stem cell, stroma score from XCELL platform; Macrophage from EPIC platform. Evidence showed that immune cells were an important influence on risk scores and tumor progression in patients with LIHC. By ssGSEA analysis, we calculated and compared the immune cell content between different risk groups and found that the expression of Macrophages (p < 0.001) and iDCs (p < 0.01) was significantly upregulated in the high-risk group compared to the low-risk group, while the expression of Mast_cells (p < 0.01) and NK_cells (p < 0.01) was significantly downregulated. Box plots of immune function indicated that Cytolytic_activity and Type_II_IFN_Reponse were significantly suppressed (p < 0.001) in the high-risk group of patients (Fig. 7B).

By scoring TME in samples from different risk groups, we found that only matrix scores differed and were higher in samples from low-risk groups (p < 0.05) (Fig. 7C). This intimated that improving the matrix score had potential benefits for patients with LIHC.

We collected 47 immune checkpoint-associated genes and performed immune checkpoint analysis on tumor samples in TCGA-LIHC. The results showed that these genes were significantly activated and upregulated (p < 0.001) in the high-risk group (Fig. 7D). Given the poorer survival outcomes in high-risk populations, immunosuppressive therapy may have a potential survival benefit for the high-risk population. We subsequently performed drug sensitivity analysis on tumor samples from the high-risk group. And the results suggested that the high-risk group was more sensitive to 32 drugs including Bortezomib, Doxorubicin, Etoposide, Embelin, and Epothilone.B (Fig. 7E). This provided a valuable reference to better guide the clinical use of medications in LIHC patients with high-risk scores.

Identification of LIHC subtypes

Applying the ConsensusClusterPlus package, we distinguished 370 tumor samples in the TCGA-LIHC cohort into 3 clusters (C1, C2, and C3), at which the clustering variable k = 3 (Fig. 8A-C). K-M plots suggested a prognostic difference between the 3 subtypes of LIHC (p < 0.001), with cluster C2 having a significantly higher OS than the other clusters (Fig. 8D). The PCA and T-SNE plots indicated that the sample distribution of different LIHC subtypes and different risk groups identified based on the whole TCGA-LIHC cohort is regionally significant (Fig. 8E-F). Sankey diagram revealed that the cluster 3 samples mostly belong to the low-risk group, while clusters 1 and 2 samples mostly belong to the high-risk group (Fig. 8G).

Immunological characterization and drug sensitivity analysis of LIHC subtypes

TME analysis suggested that the 3 clusters of LIHC do not differ in ImmuneScore, StromalScore, and ESTIMATEScore (Fig. 8H). Immune-related heat map demonstrated the degree of infiltration of immune cells from different platform sources in the 3 clusters (Fig. 8I). The results showed significant differences (p < 0.001) in the expression of some immune cells in the 3 clusters, including TIMER platform-derived Macrophage, QUANTISEQ platform-derived Monocyte, XCELL platform-derived Hematopoietic stem cell and T cell CD4 + Th2. The boxplot suggested that although most immune checkpoints had different expressions in the three clusters, the lowest activation frequently occurred in cluster C2 (Fig. 8J). It is noteworthy that samples in cluster C2 were the main source of samples from the high- and low-risk groups, so it was necessary to explore drugs that were sensitive to cluster C2. The results of drug sensitivity analysis showed that cluster C2 patients had high sensitivity to 14 drugs including AICAR, Nilotinib, Metformin, Temsirolimus, and Vorinostat (Fig. 8K).

LIHC is a worldwide health problem and deserves great attention. Chen-hao Zhang et al. calculated that there will be 906,000 newly diagnosed cases of LIHC and 83,000 LIHC-related deaths worldwide in 2022 alone, with Asian countries accounting for 72.5% and 73.3% of the mentioned data, respectively(Zhang et al. 2022a). The number of newly diagnosed cases of LIHC in developing countries, including China, has been increasing continuously. The same situation also exists in developed countries: 18,636 LIHC deaths occurred in the USA alone in 2020, and the American Cancer Society predicts this number will rise to 29,380 by 2023(Siegel et al. 2023). The global mortality rate of LIHC has continued to increase despite the fact that the manner and conditions of clinical management have improved dramatically over the past decades(Zhang et al. 2022a; Siegel et al. 2023). Notably, patients who die from LIHC are mostly in advanced stages accompanied with extensive metastases, and can only benefit from palliative care rather than surgery(Kummar and Shafi 2003). Surgical resection and liver transplantation is a privilege for LIHC patients with early-stage and can significantly prolong their survival time although the 5-year risk of recurrence remains high(Makuuchi and Sano 2004). Unfortunately, the insidious onset of LIHC with minimal symptoms, coupled with the lack of diagnostic biomarkers for early detection, often delays the diagnosis and management of LIHC to a late stage. Therefore, the development of biomarkers for early diagnosis and prognosis of LIHC would be beneficial to achieve timely and precise treatment of tumors.

The essence of LIHC progression is metastasis or in situ proliferation of tumor cells. Jamal Majidpoor and Keywan Mortezaee hold the same view: EMT and tumor cell metastasis are induced by the fracture of the tumor-stromal interface and subsequent crosstalk between tumor cells and matrix cells(Majidpoor and Mortezaee 2021). However, under normal circumstances, anoikis events will be triggered once epithelial cells separate from the ECM or adhere to abnormal ECM(Gilmore 2005). This hints that there must exist a bridge between tumor cell metastasis and the achievement of anoikis resistance. Juan Li et al. found that lncRNA dysregulation promoted tumor cell metastasis through several biological processes, including EMT, cell migration, anoikis resistance and angiogenesis(Li et al. 2016). ARlncRNAs are expected to be a novel biomarker for early diagnosis and prognosis prediction of LIHC.

In this study, information on all samples, including the tumor group (n = 374) and the control group (n = 50), were obtained from the TCGA-LIHC database. We reflected 2902 ARlncRNAs using 346 ARGs (correlation coefficient > 0.4 and p < 0.001). LIHC analysis of expression difference showed that 2143 ARlncRNAs were up-regulated and 20 ARlncRNAs were down-regulated (∣log₂FC∣>1, FDR < 0.05).

There were 2 points deserving additional attention. Firstly, we divided all tumor samples in TCGA-LIHC into a train set and test set, and the results were completely random. Second, the results obtained from the train set required additional validation, and the test and entire sets were the best choices. In this study, the preliminary results were obtained based on the training set, and the test and entire sets were used to verify the results.

After we performed the batch survival univariate Cox methods on 2163 ARlncRNAs in the train set, we screened out 27 ARlncRNAs that were strongly associated with survival prognosis (p < 0.001). LASSO regression analysis further identified an ARlncRNA signature consisting of AC100812.1, AL365295.1, AC073352.1, ELFN1-AS1, LINC00513 and MIR4435-2HG. The heatmap showed that the expression of ARlncRNAs was significantly up-regulated in tumor group. Several studies have provided partial support for our speculation and validated that upregulation of ARlncRNA expression promotes tumor metastasis. Xue Kong et al. found that AC073352.1, which is upregulated in breast cancer, is associated with the TNM progression and poor prognosis of breast cancer through promoting tumor angiogenesis and metastasis in vivo and in vitro by binding to YBX1(Kong et al. 2021). The uncontrolled upregulation of ELFN1-AS1 expression is thought to be associated with the progression of various tumors. A study by Wanguo Feng et al. suggested that ELFN1-AS1 was able to promote SBK1 expression through integration with miR-4270, ultimately leading to increased proliferation and invasion of retinoblastoma cells(Feng et al. 2021). And miR-497-3p and miR-183-3p are similar to the ELFN1-AS1 binding sponges in promoting the progression of ovarian and esophageal cancers(Zhang et al. 2020; Jie et al. 2020). In addition, ELFN1-AS1 actively participates in the metastatic activity of pancreatic, colorectal, and non-small cell lung cancers as well(Ma et al. 2021; Peng et al. 2022; Yang and Miao 2022). Chenming Zhong et al. pointed out that `MIR4435-2HG` has the potential to be a biomarker for diagnostic and prognostic use in pan-cancer(Zhong et al. 2022). Hongfei Yu et al. indicated that MIR4435-2HG is a tumor suppressor and its depletion in the stroma leads to a decrease in neutrophil content and an increase in polymorphonuclear myeloid-derived suppressor cells in tissues, which have a positive impact on the occurrence of colorectal cancer events and their progression(Yu et al. 2022). MIR4435-2HG is directly proportional to the expression level of immune cells in vivo. In the occurrence and progress of gastric cancer, macrophages are MIR4435-2HG-binding sponges which are employed to realize M2 polarization(Li et al. 2022a). Given that so many of the genomic mutations associated with cancer occur in genetic sequences that transcribe lncRNAs, we deduced that uncontrolled transcription of AC100812.1, AL365295.1, and LINC00513 probably impact tumor cell homeostasis as well(Huarte 2015).

Based on the ARlncRNAs scores, we calculated the median risk score of the prognostic model and applied it to differentiate all TCGA-LIHC samples into low- or high-risk groups. We evaluated the prognostic value of the risk model in the train, test, and entire sets separately. The results showed that the expression level of ARlncRNAs were generally upregulated in the high-risk group, and patients in the high-risk group had a worse survival status and survival prognosis than those in the low-risk group. The risk factor is more suitable for predicting the incidence of LIHC compared with other clinical features. The risk model constructed by 6 ARlncRNAs enables the maintenance of a stable predictive performance at 1, 3, and 5 years. And the prognosis of LIHC was the best in the first year.

The risk model was applicable to the prognosis prediction of patients with different clinical characteristics as well. Interestingly, the risk model was of greater prognostic value for LIHC patients with early-stage compared to those with late-stage. The results of the univariate and multivariate Cox analyses indicated that the risk score was an independent prognostic factor for LIHC. The C-index suggested that the risk model was the most accurate in assessing patient prognosis compared to other clinical traits. The risk score-based nomogram predicted the survival time of LIHC patients at 1, 3, and 5 years to 0.899, 0.801, and 0.728, respectively. The calibration and DCA curves confirmed the predictive value and performance of the nomogram. All pieces of evidence suggested that risk models constructed with 6 ARlncRNAs were essential in better achieving prognostic prediction and clinical management of LIHC.

Furthermore, we also noted that our risk model established on ARlncRNAs had a better performance in predicting the prognosis of LIHC patients compared to the prognostic models of LIHC constructed by Yutong Chen et al. and Min Deng et al. Our findings can be considered as an addition to the existing evidence in deepen our understanding of LIHC.

In our study, the survival prognosis of patients in the high-risk group was always worse than that of the low-risk group. GSEA enrichment analysis made the explanation: the KEGG signaling pathways associated with the high-risk group of LIHC were mainly involved in oocyte meiosis, pyrimidine metabolism, RNA degradation, spliceosome, and cell cycle. Zhenru Zhu et al. found that increased pyrimidine metabolism, achieved by UBE2T through stimulation of ubiquitination of Akt K63, could promote proliferation and metastasis of LIHC(Zhu et al. 2022). Jiangfeng Li et al. suggested that YTHDF2 can inhibit prostate cancer progression by triggering the degradation of some mRNAs such as LHPP and NKX3-1, and that the m6A modification site is its binding region(Li et al. 2020). Weijin Xu et al. collected 192 genes up-regulated in LIHC by meta-analysis and found that most of them were enriched in the spliceosome pathway(Xu et al. 2015). As for the cell cycle, it has been a hot spot for LIHC research. Cell cycle arrest is a central pathway for multiple drugs to exert their pharmacological effects(Abdel-Wahab et al. 2021; Niwa et al. 2022).

We performed immune cell correlation analyses on TCGA-LIHC tumor samples and synthesized the results generated from 7 platforms. The results showed that CD4 + T cells, Macrophage, and Neutrophil had the highest degree of infiltration in patients with LIHC. By ssGSEA analysis, we found that immune cells with elevated expression levels in the high-risk group were aDCs, iDCs, Macrophages, and Treg, while Macrophages were the most significant present (p < 0.001). Junfeng Wang et al. indicated that over-accumulated tumor macrophages not only originated from circulating monocytes, but were also generated from their self-renewal via tumor-derived adenosine acting in concert with autocrine granulocyte-macrophage colony-stimulating factor(Wang et al. 2021). We also observed a significant inhibition of cytolytic activity and type II IFN response in the high-risk group (p < 0.001). And this is consistent with the findings of Hideo Takahashi et al. that high depletion of cytolytic activity as an independent risk factor is associated with the worst prognosis of LIHC, leading to a decrease in the immune level and survival time of patients(Takahashi et al. 2020). TME analyses suggested that only the stromal score was statistically different between the different risk groups and that the score was lower in the high-risk group (p = 0.039). All immune checkpoints were activated and upregulated in the high-risk group (p < 0.001). Drug sensitivity analyses indicated that the high-risk group was more likely to respond to 32 drugs including Bortezomib, Doxorubicin, Etoposide, Embelin, and Epothilone. This provided guidance on the clinical use of medications in patients at high risk of LIHC.

Based on the risk score file of ARlncRNAs, we further classified TCGA-LIHC patients into 3 different subtypes by consensus cluster analysis. There were significant survival prognostic differences between the 3 subtypes of LIHC patients. PCA and t-SNE plots all suggested that the 2 risk groups and 3 LIHC subtypes we constructed were effective in differentiating between TCGA-LIHC tumor samples. Interestingly, there was no difference in TME scores between different subtypes of LIHC. The box plots indicated that LIHC subtype C2 had the lowest activation of immune checkpoint-associated genes, and that this may frustrated immune checkpoint inhibitor therapy. We remedied this deficiency by drug sensitivity analysis. The results advised that 14 kinds of drugs, including AICAR, Nilotinib, Metformin, Temsirolimus, and Vorinostat, were more effective in treating patients of the C2 subtype.

Although the data generated in this study was based entirely on in silico methods, the final results were strongly supported by the experimental data from other studies. In addition, the risk model we constructed needs to be validated with additional clinical data, such as prospective data and data from the GEO database. And the biological process regulated by ARlncRNAs will be explored in our subsequent studies.

In conclusion, we identified an anoikis-related lncRNA signature in LIHC, aiming to predict the prognosis of patients with different risk scores and to further guide their clinical dosing. ARlncRNAs revealed the immune infiltration landscape of LIHC and expands our understanding of metastasis. In an era where immune checkpoint therapy may be the only regimen that benefits LIHC patients with advanced stage, our study will contribute to the achievement of early diagnosis and personalized management of patients with LIHC.

Data Availability Statement

The datasets used and analyzed during the current study are available from the corresponding author upon reasonable request.

Author Contributions

Data collection and analysis: X.Z. and T.Z.; Writing—original draft preparation: X.Z. and T.Z.; Conception and design and final approval of the version to be published: T.Z. All authors have read and agreed to the published version of the manuscript.

Acknowledgments

This study was funded by Shanghai Municipal Health Commission [Three-year Action Plan for Developing Traditional Chinese Medicine in Shanghai, No.ZY (2018-2020)-ZYBZ-25], Jing’an District Key Clinical Specialty of Traditional Chinese Medicine Project (No. JA2020-Z003), and Jing’an District Famous Traditional Chinese Medicine Inheritance and Innovation Studio Project (No. JA2021-MLZZ005) in Shanghai.

Conflicts of Interest

The authors declare no conflict of interest.

Abdel-Wahab A-HA, Effat H, Mahrous EA, et al (2021) A Licorice Roots Extract Induces Apoptosis and Cell Cycle Arrest and Improves Metabolism via Regulating MiRNAs in Liver Cancer Cells. Nutr Cancer 73:1047–1058. https://doi.org/10.1080/01635581.2020.1783329
Chen F, Yang J, Fang M, et al (2022) Necroptosis-related lncRNA to establish novel prognostic signature and predict the immunotherapy response in breast cancer. J Clin Lab Anal 36:e24302. https://doi.org/10.1002/jcla.24302
Chen Y, Huang W, Ouyang J, et al (2023) Identification of Anoikis-Related Subgroups and Prognosis Model in Liver Hepatocellular Carcinoma. Int J Mol Sci 24:2862. https://doi.org/10.3390/ijms24032862
Collaboration GB of DLC, Akinyemiju T, Abera S, et al (2017) The Burden of Primary Liver Cancer and Underlying Etiologies From 1990 to 2015 at the Global, Regional, and National Level: Results From the Global Burden of Disease Study 2015. JAMA Oncol 3:1683–1691. https://doi.org/10.1001/jamaoncol.2017.3055
Deng M, Sun S, Zhao R, et al (2022) The pyroptosis-related gene signature predicts prognosis and indicates immune activity in hepatocellular carcinoma. Mol Med 28:16. https://doi.org/10.1186/s10020-022-00445-0
Feng W, Zhu R, Ma J, Song H (2021) LncRNA ELFN1-AS1 Promotes Retinoblastoma Growth and Invasion via Regulating miR-4270/SBK1 Axis. Cancer Manag Res 13:1067–1073. https://doi.org/10.2147/CMAR.S281536
Gilmore AP (2005) Anoikis. Cell Death & Differentiation 12:1473–1477. https://doi.org/10.1038/sj.cdd.4401723
Giuliani F, Colucci G (2009) Treatment of hepatocellular carcinoma. Oncology 77 Suppl 1:43–49. https://doi.org/10.1159/000258495
Guiu B, Garin E, Allimant C, et al (2022) TARE in Hepatocellular Carcinoma: From the Right to the Left of BCLC. Cardiovasc Intervent Radiol 45:1599–1607. https://doi.org/10.1007/s00270-022-03072-8
Han K, Kim JH (2015) Transarterial chemoembolization in hepatocellular carcinoma treatment: Barcelona clinic liver cancer staging system. World J Gastroenterol 21:10327–10335. https://doi.org/10.3748/wjg.v21.i36.10327
Hangauer MJ, Vaughn IW, McManus MT (2013) Pervasive transcription of the human genome produces thousands of previously unidentified long intergenic noncoding RNAs. PLoS Genet 9:e1003569. https://doi.org/10.1371/journal.pgen.1003569
Huarte M (2015) The emerging role of lncRNAs in cancer. Nat Med 21:1253–1261. https://doi.org/10.1038/nm.3981
Jie Y, Ye L, Chen H, et al (2020) ELFN1-AS1 accelerates cell proliferation, invasion and migration via regulating miR-497-3p/CLDN4 axis in ovarian cancer. Bioengineered 11:872–882. https://doi.org/10.1080/21655979.2020.1797281
Kong X, Li J, Li Y, et al (2021) A novel long non-coding RNA AC073352.1 promotes metastasis and angiogenesis via interacting with YBX1 in breast cancer. Cell Death Dis 12:670. https://doi.org/10.1038/s41419-021-03943-x
Kummar S, Shafi NQ (2003) Metastatic hepatocellular carcinoma. Clin Oncol (R Coll Radiol) 15:288–294. https://doi.org/10.1016/s0936-6555(03)00067-0
Lafaro KJ, Demirjian AN, Pawlik TM (2015) Epidemiology of Hepatocellular Carcinoma. Surgical Oncology Clinics of North America 24:1–17. https://doi.org/10.1016/j.soc.2014.09.001
Li C, Chen Z, Gao J, et al (2022a) MIR4435-2HG in exosomes promotes gastric carcinogenesis by inducing M2 polarization in macrophages. Front Oncol 12:1017745. https://doi.org/10.3389/fonc.2022.1017745
Li J, Meng H, Bai Y, Wang K (2016) Regulation of lncRNA and Its Role in Cancer Metastasis. Oncol Res 23:205–217. https://doi.org/10.3727/096504016X14549667334007
Li J, Xie H, Ying Y, et al (2020) YTHDF2 mediates the mRNA degradation of the tumor suppressors to induce AKT phosphorylation in N6-methyladenosine-dependent way in prostate cancer. Mol Cancer 19:152. https://doi.org/10.1186/s12943-020-01267-6
Li Z, Wang H, Deng X, et al (2022b) Construction and Validation of an Immune-Related lncRNA Prognosis Model for Thyroid Cancer. Comb Chem High Throughput Screen 25:2217–2227. https://doi.org/10.2174/1386207325666220223123917
Liu M, Zhang H, Li Y, et al (2018) HOTAIR, a long noncoding RNA, is a marker of abnormal cell cycle regulation in lung cancer. Cancer Sci 109:2717–2733. https://doi.org/10.1111/cas.13745
Llovet JM, Fuster J, Bruix J, Barcelona-Clínic Liver Cancer Group (2004) The Barcelona approach: diagnosis, staging, and treatment of hepatocellular carcinoma. Liver Transpl 10:S115-120. https://doi.org/10.1002/lt.20034
Llovet JM, Kelley RK, Villanueva A, et al (2021) Hepatocellular carcinoma. Nat Rev Dis Primers 7:6. https://doi.org/10.1038/s41572-020-00240-3
Ma G, Li G, Gou A, et al (2021) Long non-coding RNA ELFN1-AS1 in the pathogenesis of pancreatic cancer. Ann Transl Med 9:877. https://doi.org/10.21037/atm-21-2376
Majidpoor J, Mortezaee K (2021) Steps in metastasis: an updated review. Med Oncol 38:3. https://doi.org/10.1007/s12032-020-01447-w
Makuuchi M, Sano K (2004) The surgical approach to HCC: our progress and results in Japan. Liver Transpl 10:S46-52. https://doi.org/10.1002/lt.20044
Niwa AM, Semprebon SC, D’Epiro GFR, et al (2022) Salinomycin induces cell cycle arrest and apoptosis and modulates hepatic cytochrome P450 mRNA expression in HepG2/C3a cells. Toxicol Mech Methods 32:341–351. https://doi.org/10.1080/15376516.2021.2008570
Paoli P, Giannoni E, Chiarugi P (2013) Anoikis molecular pathways and its role in cancer progression. Biochimica et Biophysica Acta (BBA) - Molecular Cell Research 1833:3481–3498. https://doi.org/10.1016/j.bbamcr.2013.06.026
Peng F, Wang R, Zhang Y, et al (2017) Differential expression analysis at the individual level reveals a lncRNA prognostic signature for lung adenocarcinoma. Mol Cancer 16:98. https://doi.org/10.1186/s12943-017-0666-z
Peng S, Luo Y, Chen L, et al (2022) lncRNA ELFN1-AS1 enhances the progression of colon cancer by targeting miR-4270 to upregulate AURKB. Open Med (Wars) 17:1999–2012. https://doi.org/10.1515/med-2022-0582
Shi T, Gao G, Cao Y (2016) Long Noncoding RNAs as Novel Biomarkers Have a Promising Future in Cancer Diagnostics. Dis Markers 2016:9085195. https://doi.org/10.1155/2016/9085195
Siegel RL, Miller KD, Wagle NS, Jemal A (2023) Cancer statistics, 2023. CA Cancer J Clin 73:17–48. https://doi.org/10.3322/caac.21763
Song J, Liu Y, Liu F, et al (2021) The 14-3-3σ protein promotes HCC anoikis resistance by inhibiting EGFR degradation and thereby activating the EGFR-dependent ERK1/2 signaling pathway. Theranostics 11:996–1015. https://doi.org/10.7150/thno.51646
Taddei ML, Giannoni E, Fiaschi T, Chiarugi P (2012) Anoikis: an emerging hallmark in health and diseases. J Pathol 226:380–393. https://doi.org/10.1002/path.3000
Takahashi H, Kawaguchi T, Yan L, et al (2020) Immune Cytolytic Activity for Comprehensive Understanding of Immune Landscape in Hepatocellular Carcinoma. Cancers (Basel) 12:1221. https://doi.org/10.3390/cancers12051221
Vogel A, Meyer T, Sapisochin G, et al (2022) Hepatocellular carcinoma. The Lancet 400:1345–1362. https://doi.org/10.1016/S0140-6736(22)01200-4
Wang J, Wang Y, Chu Y, et al (2021) Tumor-derived adenosine promotes macrophage proliferation in human hepatocellular carcinoma. J Hepatol 74:627–637. https://doi.org/10.1016/j.jhep.2020.10.021
Xia S, Wu J, Zhou W, et al (2021) HRC promotes anoikis resistance and metastasis by suppressing endoplasmic reticulum stress in hepatocellular carcinoma. Int J Med Sci 18:3112–3124. https://doi.org/10.7150/ijms.60610
Xu W, Huang H, Yu L, Cao L (2015) Meta-analysis of gene expression profiles indicates genes in spliceosome pathway are up-regulated in hepatocellular carcinoma (HCC). Med Oncol 32:96. https://doi.org/10.1007/s12032-014-0425-6
Xu Y, Zhang X, Hu X, et al (2018) The effects of lncRNA MALAT1 on proliferation, invasion and migration in colorectal cancer through regulating SOX9. Mol Med 24:52. https://doi.org/10.1186/s10020-018-0050-5
Yang B, Miao S (2022) lncRNA ELFN1-AS1 predicts poor prognosis and promotes tumor progression of non-small cell lung cancer by sponging miR-497. Cancer Biomark 34:637–646. https://doi.org/10.3233/CBM-210393
Yu H, Chen C, Han F, et al (2022) Long Noncoding RNA MIR4435-2HG Suppresses Colorectal Cancer Initiation and Progression By Reprogramming Neutrophils. Cancer Immunol Res 10:1095–1110. https://doi.org/10.1158/2326-6066.CIR-21-1011
Zhang C, Lian H, Xie L, et al (2020) LncRNA ELFN1-AS1 promotes esophageal cancer progression by up-regulating GFPT1 via sponging miR-183-3p. Biol Chem 401:1053–1061. https://doi.org/10.1515/hsz-2019-0430
Zhang C-H, Cheng Y, Zhang S, et al (2022a) Changing epidemiology of hepatocellular carcinoma in Asia. Liver Int 42:2029–2041. https://doi.org/10.1111/liv.15251
Zhang P, Li H, Cheng X, Wang W (2022b) Comprehensive Analysis of Immune Cell Infiltration of m6a-Related lncRNA in Lung Squamous Cell Carcinoma and Construction of Relevant Prognostic Models. Biomed Res Int 2022:9139823. https://doi.org/10.1155/2022/9139823
Zhong C, Xie Z, Zeng L-H, et al (2022) MIR4435-2HG Is a Potential Pan-Cancer Biomarker for Diagnosis and Prognosis. Front Immunol 13:855078. https://doi.org/10.3389/fimmu.2022.855078
Zhong X, Rescorla FJ (2012) Cell surface adhesion molecules and adhesion-initiated signaling: Understanding of anoikis resistance mechanisms and therapeutic opportunities. Cellular Signalling 24:393–401. https://doi.org/10.1016/j.cellsig.2011.10.005
Zhu Z, Cao C, Zhang D, et al (2022) UBE2T-mediated Akt ubiquitination and Akt/β-catenin activation promotes hepatocellular carcinoma development by increasing pyrimidine metabolism. Cell Death Dis 13:154. https://doi.org/10.1038/s41419-022-04596-0

No competing interests reported.

Download PDF

Version 1

posted

You are reading this latest preprint version

An Anoikis-related LncRNA Signature: Integrated Analysis of Immune Infiltration Landscape in Patients with Liver Hepatocellular Carcinoma and Immunotherapy

Status:

Version 1

Abstract

Figures

Introduction

Materials and Methods

Results

Discussion

Declarations

References

Additional Declarations

Status:

Version 1