Transcriptomic Profile of Lymphovascular Invasion, a Known Risk Factor of Pancreatic Ductal Adenocarcinoma Metastasis

Lymphovascular invasion (LVI) is an aggressive pathologic feature and considered a risk factor for distant metastasis. We hypothesized that pancreatic ductal adenocarcinomas (PDACs) with LVI are associated with shorter survival, as well as aggressive cancer biology and lymphangiogenesis in transcriptomic analysis. Utilizing the cancer genome atlas (TCGA)-PDAC cohort, we found that positive LVI was significantly associated with positive perineural invasion (PNI) (p = 0.023), and higher American Joint Committee on Cancer (AJCC) T (p = 0.017) and N (p < 0.001) categories. Furthermore, positive LVI was associated with shorter overall survival (OS) (p = 0.014) and was an independent risk factor of poor OS. Although there was no association between LVI status and lymphangiogenesis, epithelial-mesenchymal transition (EMT), or metastasis-related genes, Gene Set Enrichment Analysis revealed a strong association with cell-proliferation-related gene sets such as mitotic spindles (Normalized enrichment score (NES) = 1.76, p = 0.016) and G2/M checkpoints (NES = 1.75, p = 0.036), as well as with transforming growth factor beta (TGF-beta) signaling (NES = 1.61, p = 0.043), which is a known mechanism of metastasis in PDACs. In conclusion, positive LVI was an independent risk factor of poor OS in PDACs. We found that PDACs with LVI were possibly associated with accelerated cell proliferation and enhanced TGF-beta signaling independent of lymphangiogenesis. Transcriptomic profiling elucidates more precise tumor biology of LVI-positive PDACs.


Introduction
Pancreatic ductal adenocarcinoma (PDAC) will become the second leading cause of cancer-related deaths in the United States by 2030 [1,2]. The exact reason for the increased incidence of PDACs is

LVI Status Was Associated with Worse Prognosis in PDAC
Survival analysis revealed that, although it did not reach statistical significance, patients with positive LVI tumors tended to have a shorter median disease-free survival (DFS) compared with patients without LVI (median DFS time: 12.2 vs. 20.4 months, p = 0.069, Figure 2A). Furthermore, patients with LVI-positive tumors had a significantly shorter median overall survival (OS) compared with patients with LVI-negative tumors (median OS time: 17

LVI Status Was Associated with Worse Prognosis in PDAC
Survival analysis revealed that, although it did not reach statistical significance, patients with positive LVI tumors tended to have a shorter median disease-free survival (DFS) compared with patients without LVI (median DFS time: 12.2 vs. 20.4 months, p = 0.069, Figure 2A). Furthermore, patients with LVI-positive tumors had a significantly shorter median overall survival (OS) compared with patients with LVI-negative tumors (median OS time: 17.0 vs. 22.5 months, p = 0.014, Figure 2B). Univariate analysis showed that residual tumor status (p = 0.021), adjuvant radiation (p = 0.002), and positive LVI status (p = 0.011) were associated with OS. Multivariate analysis revealed that positive LVI was independently associated with shorter OS (p = 0.012, Table 2). Positive LVI was found to affect clinical outcome, which is in agreement with the previous reports [14,15]. LVI was significantly associated with shorter OS in the TCGA cohort, demonstrating that the TCGA-PDAC cohort is a reliable cohort that follows a similar trend as the previous reports in the literature. 0.002), and positive LVI status (p = 0.011) were associated with OS. Multivariate analysis revealed that positive LVI was independently associated with shorter OS (p = 0.012, Table 2). Positive LVI was found to affect clinical outcome, which is in agreement with the previous reports [14,15]. LVI was significantly associated with shorter OS in the TCGA cohort, demonstrating that the TCGA-PDAC cohort is a reliable cohort that follows a similar trend as the previous reports in the literature.    Given the previous reports that demonstrated higher expression of vascular and lymphatic angiogenesis-related genes in the LVI-positive tumors in breast and colon cancer [30,34,35], we speculated that vascular and lymphatic angiogenic genes might be upregulated in the LVI-positive PDACs. However, contrary to the previous reports, LVI was not associated with the expressions of lymphangiogenesis-related genes, such as PDPN, neuropilin-2 (NRP2), vascular endothelial growth factor (VEGF) C, Angiopoietins (ANGPT1 and ANGPT2), hypoxia-inducible factor (HIF)-1α, platelet-derived growth factor (PDGF)-BB, fibroblast growth factor (FGF), lymphatic vessel endothelial hyaluronan receptor (LYVE)-1, and prospero homeobox protein (PROX) 1 ( Figure 3A), as well as vascular-angiogenesis-related genes, such as VEGFA and VEGFB ( Figure 3B). Further, with gene set enrichment analysis (GSEA), the angiogenesis gene set was not enriched in the LVI-positive PDACs (normalized enrichment score (NES) = 1.17, p = 0.160, Figure 3C). Next, we compared the expression levels of genes involved in ECM degradation, such as matrix metalloproteinase (MMP)1, MMP9, and MMP14, which were previously reported to be elevated in LVI-positive breast cancers [31][32][33]. However, none of the MMP expressions were significantly different between the LVI-positive and negative PDACs ( Figure 3D). enrichment analysis (GSEA), the angiogenesis gene set was not enriched in the LVI-positive PDACs (normalized enrichment score (NES) = 1.17, p = 0.160, Figure 3C). Next, we compared the expression levels of genes involved in ECM degradation, such as matrix metalloproteinase (MMP)1, MMP9, and MMP14, which were previously reported to be elevated in LVI-positive breast cancers [31][32][33]. However, none of the MMP expressions were significantly different between the LVI-positive and negative PDACs ( Figure 3D).

Positive LVI Was Associated with Promoted Cell Cycles and Enhanced Transforming Growth Factor (TGF)-Beta Signaling in PDAC
Given that the positive LVI status was associated with increased lymph node metastases and a short median OS, we hypothesized that there are associations between LVI and aggressive cancer features, such as epithelial-mesenchymal transition (EMT), tumor inflammation, metabolism, cell proliferation, cancer stem cell maintenance, and several other important signaling pathways in PDACs. As Jones et al. reported that there were several core signaling pathways altered in most PDACs, including transforming growth factor beta (TGF-beta), Notch signaling, WNT-beta catenin, and Hedgehog pathways [36], we suspected that a positive LVI status might associate with these pathways. GSEA revealed that cell-cycle-related gene sets, including mitotic spindle (NES = 1.76, p = 0.016, Figure 4A) and G2/M checkpoint (NES = 1.75, p = 0.036, Figure 4B), were significantly enriched in the LVI-positive PDACs, whereas the LVI-positive PDACs did not enrich EMT (NES = 1.19, p = 0.344, Figure 4C), apoptosis (NES = 1.03, p = 0.404, Figure 4D), glycolysis (NES = 1.22, p = 0.237, Figure 4E), or the inflammatory response gene set (NES = 0.60, p = 0.901, Figure 4F). Additionally, there were no associations with cancer stem cell regulation or maintenance genes, including Forkhead box protein Q1 (FOXQ1), B cell-specific Moloney murine leukemia virus integration site 1 (BMI1), and Enhancer of zeste homolog 2 (EZH2) [37] (Supplemental Figure S1). Among the reported altered core signaling pathways, the TGF-beta signaling gene set was enriched in the tumors with positive LVI (NES = 1.61, p = 0.043, Figure 5A). However, Notch signaling (NES = 0.97, p = 0.498, Figure 5B), WNT-beta catenin signaling  Figure 5C), and Hedgehog signaling (NES = 0.80, p = 0.718, Figure 5D) were not enriched in the LVI-positive tumors. We further investigated neuropilin-1 (NRP1) expression, as it acts as a co-receptor for TGF-beta. However, the NRP1 expression was not associated with LVI status ( Figure 5E). This may be because NRP1 interacts not only with TGF-beta, but also with multiple other growth factors [38]. These findings suggest that positive LVI is associated with promoted cell cycles and enhanced TGF-beta signaling pathways, possibly leading to a worse prognosis with higher tumor growth and potential to metastasize to lymph nodes. 0.237, Figure 4E), or the inflammatory response gene set (NES = 0.60, p = 0.901, Figure 4F). Additionally, there were no associations with cancer stem cell regulation or maintenance genes, including Forkhead box protein Q1 (FOXQ1), B cell-specific Moloney murine leukemia virus integration site 1 (BMI1), and Enhancer of zeste homolog 2 (EZH2) [37] (Supplementary Figure S1). Among the reported altered core signaling pathways, the TGF-beta signaling gene set was enriched in the tumors with positive LVI (NES = 1.61, p = 0.043, Figure 5A). However, Notch signaling (NES = 0.97, p = 0.498, Figure 5B), WNT-beta catenin signaling (NES = 0.78, p = 0.764, Figure 5C), and Hedgehog signaling (NES = 0.80, p = 0.718, Figure 5D) were not enriched in the LVI-positive tumors. We further investigated neuropilin-1 (NRP1) expression, as it acts as a co-receptor for TGF-beta. However, the NRP1 expression was not associated with LVI status ( Figure 5E). This may be because NRP1 interacts not only with TGF-beta, but also with multiple other growth factors [38]. These findings suggest that positive LVI is associated with promoted cell cycles and enhanced TGF-beta signaling pathways, possibly leading to a worse prognosis with higher tumor growth and potential to metastasize to lymph nodes.

Discussion
To our knowledge, this is the first study evaluating the association between pathologically proven LVIs and transcriptomic information in a PDAC cohort. Some of the pathological information, such as LVI or PNI, is not included in clinicopathological data in the TCGA cohort. Howev-

Discussion
To our knowledge, this is the first study evaluating the association between pathologically proven LVIs and transcriptomic information in a PDAC cohort. Some of the pathological information, such as LVI or PNI, is not included in clinicopathological data in the TCGA cohort. However, with the Text Information Extraction System (TIES), we were able to access the de-identified original pathology reports on each TCGA patient and link pathological information to the transcriptomic data manually. This methodology enabled us to complete this study.
Similar to our previous study on the association between enhanced vascularity in PDACs and better survival, we investigated the association between LVI-positive PDACs and lymphangiogenesis and ECM degradation, utilizing transcriptomic profiling analysis. We found that PNI status and AJCC T and N categories were associated with LVI status in the PDACs. Using the TCGA cohort, we confirmed that positive LVI status was significantly associated with poor survival in PDAC and an independent risk factor for OS, echoing the previous reports [14,15]. Further, transcriptomic investigation revealed that the PDACs with positive LVI were associated with enhanced cell cycles and promoted TGF-beta signaling pathways, but not with EMT or lymphangiogenesis.
We found that that positive LVI tumors were associated with a higher AJCC T category and increased lymph node metastases, which is in agreement with previous studies [5,7,14]. Epstein et al. reported that PDACs with positive LVI have higher potential to metastasize to the regional lymph nodes through the lymphatic channels, yet LVI is also a significant risk factor for survival independent of lymph node metastasis [14]. Although the tumor size between the LVI-positive and negative tumors was similar, LVI-positive tumors appeared to be more proliferative. As commonly seen in clinical practice, small PDACs do have the ability to metastasize to lymph nodes as well as distant organs, which makes it conceivable that tumor size is only one aspect of tumor aggressiveness. LVI is a complicated invasion process, requiring tumor cells to detach themselves from the growing primary mass, followed by invasion and migration through the ECM towards vascular walls. In PDACs, VEGFs were reported to associate with enhanced angiogenesis and lymphangiogenesis [39,40]. Additionally, in breast cancers, a strong correlation has been reported between LVI and genes controlling ECM degradation, angiogenesis, neovascularization, and tumor necrosis factor [28,29,31]. Due to scarce evidence of transcriptomic characterization with LVI in PDACs, we investigated the association between LVI status and these previously reported genes, which unexpectedly revealed no association in PDACs. Furthermore, EMT is known to be one of the major mechanisms of PDAC metastasis, as the transformation is associated with increased mobility of the cancer cells [41]. Transcriptomic profiles in the present study demonstrated no association between EMT and LVI status; however, this may be because EMT activation might have occurred in an early phase during pancreatic cancer initiation [42], whereas LVI-positive tumors are in advanced PDACs. Positive LVI was associated with enhanced cell cycles, such as mitotic spindle and G2/M checkpoint gene upregulation. These results suggest that cancer cells with enhanced cell cycles and rapid proliferation [43] lead to lymphovascular invasion and subsequent lymph node metastases. Interestingly, our group recently found that breast cancer with LVI similarly demonstrated enhanced cellular proliferation independent of lymphangiogenesis, which suggests that LVI might be simply a marker of tumor proliferation with worse tumor biology, rather than increased lymphangiogenesis [44]. We cannot help but speculate that EMT and lymphangiogenesis may not be as clinically relevant as tumor proliferation based on these results.
TGF-beta signaling was associated with positive LVI status in the present study. TGF-beta signaling is one of several core signaling pathways commonly altered in PDACs [36], and enhanced expression is associated with poor prognoses [45,46]. TGF-beta mediates interaction between pancreatic stellate cells and cancer cells, producing a rich stromal reaction, which is one of the major characteristics of PDACs [47]. TGF-beta signaling also regulates immune cells in TME, inhibits the functions of natural killer (NK) and CD8+ T cells, induces pro-oncogenic regulatory T cells (Treg), and inhibits B cell proliferation [48]. Additionally, Bierie highlighted a paracrine mechanism of TGF-beta, interacting between various cell populations in TME and playing a significant role in promoting cell proliferation by secreting growth factors and increasing metastatic potential [49]. Interestingly, TGF-beta is also known to have a dual role in PDAC progression: a pro-oncogenic role and an anti-oncogenic role [45]. While TGF-beta acts as a tumor suppressor in early-stage pancreatic cancer by inhibiting epithelial cell cycle progression and promoting apoptosis, it becomes pro-oncogenic in later stages by promoting EMT, inducing pro-oncogenic immune cells in TME, increasing cell motility, and increasing distant metastasis [45,50].
These novel findings should be noted in the context of the limitations of this study. First, we utilized solely publicly available TCGA datasets. We are aware that our findings and conclusions would have been stronger with a validation cohort. Unfortunately, given the low incidence of pancreatic cancer despite its high mortality, a cohort including transcriptomic and clinicopathological data was not available, despite our extensive search. Second, there was some missing information from the original reports in TIES. However, we were able to retrieve approximately 85% of LVI and PNI statuses. Further, pathological information for LVI and PNI was obtained from the original reports through TIES, and these reports were not validated or standardized [51]. Third, whereas TCGA has significant benefits, with a large amount of transcriptomic information associated with patient survival, it also has limitations, such as a relatively short follow-up duration, as well as no detailed information available for neoadjuvant or adjuvant chemotherapy. Additionally, the cohort only contains surgically resected specimens, not metastatic sites, which could possibly make our conclusion difficult to generalize. In addition, we did not discuss developing definitive biomarkers to differentiate patients with a high risk of LVI, since the aim of this study was to elucidate the underlying cancer biology of LVI. Development of a new biomarker to differentiate PDAC patients with worse outcomes will be interesting and will be pursued in a different manuscript. Lastly, this study does not include any in vitro or in vivo experiments; therefore, all of our findings were observational associations. To this end, our study did not prove any causal mechanisms; however, since this study elucidated the underlying biology of LVIs using transcriptomic profiles from actual human pancreatic cancer, this study is more precise than in vivo or in vitro experiments in terms of reflecting human cancer. Despite these limitations, we believe that this study still provided valuable underlying aggressive biological information on LVI-positive PDACs.

Data Acquisition from TCGA
Clinicopathological data and gene expression data from RNA sequences were obtained from the TCGA pancreatic cancer cohort (PAAD) through cBioPortal [52,53], which contained 185 primary pancreas tumors, as described in the past [54][55][56]. One hundred fifty-four patients were registered as "Pancreas-Adenocarcinoma Ductal Type" within OS data [57,58]. The median follow-up period was 15 months (Inter-quartile range (IQR): 8-21 months). Gene expression data were utilized after being log2 transformed. LVI and PNI information associated with the TCGA cohort was manually obtained via the TIES system that attaches the pathological reports to the TCGA cohort (http://ties.dbmi.pitt.edu/#) through the Roswell Park Comprehensive Cancer Center. Given that TCGA is a publicly available, de-identified database, the Institutional Review Board (IRB) was exempted.

Survival Analysis
DFS was defined as the time from the date of surgery to the date of disease progression. OS was defined as the time from the date of surgery to the date of death by any cause. In order to obtain hazard ratios (HRs) and 95% confidence intervals (CIs), univariate and multivariate analyses for OS were conducted, using age, gender, tumor size, AJCC Staging T and N categories, histological grade, residual tumor status, history of adjuvant radiation treatment, PNI, and LVI. AJCC pathological stage IA and IB and stage IIA and IIB were simplified to stage I and stage II, respectively. Some of the parameters were dichotomized as follows: age below 65 and age 65 and above, tumor size below 3.5 cm and greater than or equal to 3.5 cm (the median value), AJCC T category T1 + T2 and T3 + T4, pathological AJCC stage I + II and III + IV, residual tumor status R1 + R2 and R0, and histologic grade G3 (poorly differentiated) and G1 + G2 (well-differentiated + moderately-differentiated). Multivariate analysis was performed with variables with a p-value < 0.05 on univariate analyses.

Statistical Analysis
Statistical comparisons of the clinicopathological parameters were performed by a Chi square test. Continuous values were compared by a Student's t-test. All statistical analyses were performed using JMP 14.2 (SAS, Cary, NC, USA) and R software (http://www.r-project.org/) together with Bioconductor (http://bioconductor.org/). A p-value < 0.05 (two-tailed) was considered to be statistically significant.

Conclusions
In conclusion, a positive LVI status was associated with shorter OS in PDACs, possibly due to accelerated cell proliferation and enhanced TGF-beta signaling, leading to faster tumor growth and higher lymph node metastasis potential. Transcriptomic profiling from human pancreatic cancer specimens elucidates more precise tumor biology in LVI-positive PDACs. Funding: This work was supported by National Institute of Health (NIH) grant R01CA160688 to K.T., and National Cancer Institute (NCI) grant P30CA016056 and U24CA232979 involving the use of the Roswell Park Comprehensive Cancer Center's Bioinformatics and Biostatistics Shared Resources. Additionally, this research used the TIES, which is supported by NCI grant U24 CA180921.

Conflicts of Interest:
The authors declare no conflict of interest.