The public and de-identified data of primary solid tumors from TCGA database by cBioPortal for Cancer Genomics (https://www.cbioportal.org, accessed in March 25, 2020), including genomic, demographic, clinicopathological and prognostic data were retrieved (14). The genomic data consisted of mRNA-Seq expression data, which was generated using Illumina HiSeq V2 platform. The mRNA-Seq data were processed and normalized using RNA-Seq by expectation maximization, and transformed into log 2 values for analysis (15). Our pan-cancer analysis only focused on Mu Opioid Receptor mRNA expression in 18 common primary solid tumors ( Table S1 ) and any patients whose MOR mRNA expression was not available were excluded.

The MOR mRNA-Seq expression data (symbol: OPRM1; gene ID: 4988) of 18 common solid cancer types were retrieved from the TCGA database. Initially, we divided all patients into positive versus negative subgroups by median MOR mRNA expression values of each cancer type (16). Then we found MOR mRNA expression was at a low level, and most median values were zero. Finally, we defined patients with zero MOR expression as negative. Demographic data included age, gender and race. Age was classified as young (under 60 years old), old (over 60 years old) and unknown subgroups (the median of age was 60 years old). Gender was classified as male, female and unknown. Race was classified as Caucasian, African and others. Clinicopathological data were the American Joint Committee on Cancer (AJCC) Tumor-node-metastasis (TNM) stages and histological grade. AJCC pathologic TNM stage was classified as stage I, stage II, stage III, stage IV and others. AJCC pathologic T stage was classified as T1, T2, T3, T4 and others. N stage was classified as N0, N1, N2, N3 and others. M stage was classified as M0, M1 and others. Histological grade was classified as low grade, high grade and others. In this study, when investigating the prognostic features of MOR, we focused on overall survival (OS), progression-free survival (PFS), disease-free survival (DFS) and disease-specific survival (DSS).

We also analyzed the association between MOR mRNA expression and tumor microenvironmental immune types (TMIT) and assessed the immunogenic features of the MOR. According to previous studies (17), the TMIT classification was divided into high or low expressions based on the expression of PD-L1 and CD8A. The type I represents high PD-L1 and high CD8A expression, type II with low PD-L1 and low CD8A expression, type III as high PD-L1 and low CD8A expression, and type IV of low PD-L1 and high CD8A expression. To explain, TMIT type II implies decreased infiltration of CD8+ T cells in the tumor microenvironment and decreased expression of PD-L1 in cancer cells, which represents a poor response to PD-L1 therapy.

Continuous variables with normal distribution were described as mean and standard deviation. Continuous variables with skewed distribution were described as median, first quartile and third quartile. Categorical variables were described as frequency and percentage. Pearson’s Chi-squared test or Fisher’s exact test was used to detect the statistical significance for categorical variables of demographic and clinicopathological features between MOR subgroups, as well as the association between TMIT and MOR mRNA expression. Independent Student’s t-test was used to detect a statistical significance for continuous variables with normal distribution and homogeneity of variances between MOR subgroups; Otherwise, Kruskal-Wallis test was used. Binomial logistic regression models were used to detect associations between MOR mRNA expression and binary clinicopathological features. Multinomial logistic regression models were used to detect the associations between MOR mRNA expression and polyfactorial clinicopathological features. The greater odds ratio (OR) value indicates the more advanced stage and grade of cancer. The log-rank test and Kaplan-Meier estimator were used to screen significant prognostic factors that were associated with survival outcomes. Glioblastoma multiforme (GBM) and skin cutaneous melanoma (SKCM) were excluded from the disease-free survival analysis, because no disease-free survival data was available for GBM and SKCM. After adjustment of Cox regression model, the hazard ratio (HR) for each prognostic factor was calculated. The greater HR value suggested the greater risk of death. To eliminate potential disequilibrium caused by confounding factors, subgroup analysis and propensity score matching (PSM) were done. All statistical analysis was done by R statistical software (version 3.6.2, released on February 29, 2020). A Two-sided P value < 0.05 was considered to be of statistical significance.

There were 7,274 patients with 18 common solid cancer types included into this study; Of those, 6126 patients were with negative MOR mRNA expression and 1,112 patients with positive MOR mRNA expression ( Table 1 ). There were 3,630 young patients with age up to 60 years old, 3,612 male patients, and 5, 222 Caucasian patients ( Table 1 ). The majority of patients belonged to stage I (1755), N0 stage (3247) and M0 stage (3872) ( Table 1 ). Generally, the positive rates of MOR mRNA expression varied from 2.7% to 50.2% across different cancer types ( Table S1 ). The clinicopathological data were missing in some cancer types, and these cancer types were excluded in subsequent analysis. Compared with normal tissues, MOR mRNA was generally expressed in solid cancers, but the mRNA expression level was relatively low ( Figure S1 ). The colon adenocarcinoma (COAD) and rectal adenocarcinoma (READ) had similar embryological and histological features and they were classified as one group (colon and rectal adenocarcinoma, COREAD) in the following analysis.

We explored the association between clinicopathological features (including AJCC TNM stage, histological grade, T stage, N stage and M stage) and mRNA expression of MOR across 18 major solid cancer types. Patients with positive MOR mRNA expression tended to have more advanced AJCC TNM stage, T stage, N stage and M stage when compared with patients with negative MOR mRNA expression ( Figure 1 ). There was no visible relationship between MOR expression and histological grade. In the following analysis, we combined N1, N2 and N3 together.

After adjustment of multivariate logistic regression analysis, as shown in Table 2 , we found that the association between T stage and MOR mRNA expression remained statistically significant (adjusted OR, 1.176; 95% confidence interval [CI], 1.022-1.354; P=0.024), as well as the association between M stage and MOR mRNA expression (adjusted OR, 1.548; 95% CI, 1.095-2.189; P=0.013), but association between N stage and MOR was not significant statistically (adjusted OR, 1.145; 95% CI, 0.975-1.346; P=0.101). In the sensitivity analysis, we removed one cancer type at each time and then re-analyze the association between clinicopathological features and MOR mRNA expression. We found that the association between T stage, M stage and MOR mRNA expression was stable ( Figure S2 ). Subgroup analysis showed that there was an association trend between clinicopathological features and MOR mRNA expression (most adjusted OR>1.0), but the statistical difference is not significant (most P>0.05), ( Figure S3 ).

There were 43 patients (0.6%) with missing OS data, 43 patients (0.6%) with missing PFS data, 3529 patients (48.5%) with missing DFS data, and 275 patients (3.8%) with missing DSS data. These patients were excluded from the corresponding survival analysis. The median follow-up periods were 23.87 months (13.12-46.45 months) for OS, 19.79 months (10.00-39.45 months) for PFS, 23.74 months (13.71-43.73 months) for DFS, and 23.87 months (13.25-46.32 months) for DSS.

The prognosis of patients with positive MOR mRNA expression was worse than patients with negative MOR expression (P<0.001 for OS, PFS, DFS and DSS) ( Figure 2 ). The poor prognosis of positive MOR mRNA expression was observed in most cancer types ( Figures S4-S7 ). Multivariate Cox regression models identified positive MOR mRNA expression as a significant prognostic factor in all cancer types (HR 1.347, 95% CI 1.200-1.512, P<0.001 for OS; HR 1.359, 95% CI 1.220-1.513, P<0.001 for PFS; HR 1.269, 95% CI 1.016-1.585, P<0.001 for DFS; HR 1.474, 95% CI 1.284-1.693, P<0.001 for DSS) ( Figure S8 ). Sensitivity analysis confirmed the prognostic significance of MOR mRNA expression ( Figure S8 ). After removing liver hepatocellular carcinoma (LIHC) from the study cohort, the effect of positive MOR mRNA expression on DFS turned more remarkable (HR 1.539, 95% CI 1.218−1.944, P<0.001) ( Figure S8 ). We applied PSM to generate balanced data from the raw data, in which the demographics differences were eliminated between the two groups ( Table S2 ). The adverse effect of positive MOR expression on prognosis still existed (P<0.001 for OS, PFS, DFS and DSS) with the balanced data ( Figure S9 ).

We grouped all the patients into four TMIT groups according to the expression levels of CD8A and PD-L1. Among the included patients, 35.0% were classified as TMIT I, with high PD-L1 expression and high CD8A expression. High PD-L1 expression indicated a favorable response to the PD-L1 therapies. High CD8A expression represented a high proportion of CD8+ CTLs in the tumor microenvironment. The proportions of TMIT II (low PD-L1/low CD8A), III (high PD-L1/low CD8A), and IV (low PD-L1/high CD8A) were 35.0%, 15.0% and 15.0%, respectively ( Table S3 ).

We analyzed the relationship between MOR mRNA expression and TMIT. As patients with positive MOR mRNA expression were likely to be classified as TMIT II, while patients with negative MOR expression were likely to be classified as TMIT I ( Figure 3A ). Patients of TMIT II had higher MOR mRNA expression, while patients of TMIT I had lower MOR expression ( Figures 3B, D and Figure S10 ). This pattern remained in breast invasive carcinoma (BLCA), cervical squamous cell carcinoma (CESC), head and neck squamous cell carcinoma (HNSC), lung squamous cell carcinoma (LUSC). The expression levels of PD-L1 and CD8A were positively correlated overall and in most cancer types ( Figure 3C and Figure S11 ).

To summarize, MOR mRNA overexpression in solid tumors was associated with advanced cancer and poor survival, which may be related to the tumor microenvironment ( Figure 4 ). In the tumor microenvironment, increased MOR expression indicated decreased PD-L1 expression on cancer cells and decreased CD8+ T cell infiltrations, demonstrating poor response to PD-L1 therapy ( Figure 4 ).