IntroductionMultiple myeloma (MM) is a hematologic malignancy characterized by the abnormal proliferation and accumulation of plasma cells (PCs) in the bone marrow (BM) [1]. MM cells present genetic and molecular aberrations, which are defined as primary and secondary genetic events [2]. Among the secondary cytogenetic alterations, copy number alteration (CNA) of the 1q21 region (1q21+) is one of the most frequent events [3, 4]. The incidence and copy number phenotype of the 1q21 region increase during the progression and relapse of the disease. In fact, additional copies of 1q21 can be detected in around 40% of smoldering (SMM) and newly diagnosed MM (NDMM) and around 70% of relapse-refractory MM (RRMM) [5, 6].As a gold standard, 1q21+ can be detected by fluorescent in situ hybridization (FISH). Gain (three copies) or amplification [amp (four or more copies)] of the 1q21 region occurs as a progression event and provides a fitness advantage to a particular sub-clone over other populations, contributing to tumor progression and relapse [7]. Particularly, MM patients with four or more copies of 1q21 are associated with high-risk features and poor response to standard therapies [4, 8]. Hanamura et al. showed that patients with 1q21 amp at diagnosis tended to have a more aggressive clinical course than those with 1q21 gain [5]. Moreover, the presence of 1q21+ has also been correlated with a reduced response and/or resistance to proteasome inhibitors (PIs) treatment, including bortezomib (BTZ) and carfilzomib (CFZ) [8,9,10] in MM patients.The 1q21 amplicon is known to contain several genes that have been suggested to drive disease aggressiveness in 1q21+ patients [11,12,13,14,15]. However, the exact genes driving the 1q21+ have not been fully characterized. Furthermore, the mechanisms underlying the high-risk disease progression associated with 1q21+ have not been identified. Nevertheless, the identification in the 1q21 region of druggable targets is an emerging unmet medical need in MM patients for a personalized approach.For this reason, in our study, we initially analyzed the gene expression of CD138+ PCs obtained from MM BM aspirates. This analysis allowed us to identify a set of genes within the 1q21 region, with copy number-driven expression. Among the genes most expressed in patients with 1q21 CNA compared to controls, we identified and focused on PYGOPUS2 (PYGO2), a transcriptional factor involved in the Wnt/β-catenin signaling pathway. Across multiple solid malignancies, including glioma, prostate, liver, lung and breast cancer, PYGO2 is frequently overexpressed and correlates with enhanced proliferation, epithelial-mesenchymal transition and adverse clinical outcome [16,17,18]. Particularly, in prostate cancer, PYGO2 has been found to drive disease progression and metastasis [19].The PYGO2 protein, composed of 406 amino acids with a molecular mass of 41 kDa, has two conserved domains, an N-terminal homology domain (NHD) and a C-terminal plant homology domain (PHD) zinc finger motif. Studies carried out in central nervous system tumors showed that deletion of the NHD domain is associated with 50% reduction of transcriptional activity [20]. The PHD domain has an important role in transcriptional activation, since PYGO2 is implicated in chromatin remodeling and in binding to methylated residues on lysine 4 of histone H3 (H3K4me3) [21]. The plant PHD finger is responsible for the binding, through adaptor proteins, to the N-terminal domain of β-catenin [21].Particularly, PYGO2 is a nuclear protein involved in the TCF/β-catenin transcriptional complex, and it was first identified as a significant component of the canonical Wnt/β-catenin pathway, forming a linear chain with β-catenin via BCL9. This chain anchors β-catenin in the nucleus and activates transcription of target genes in the Wnt/β-catenin pathway [22, 23]. The target genes regulated by PYGO2 are Cyclin D1 (CCND1), C-myc, Cyclin A, CD44 and multidrug resistance 1 (MDR1/ABCB1), which all play important roles in tumorigenesis, tumor progression, and prognosis. This emphasizes the implication and relevance of the PYGO2 protein in cell growth and proliferation [18, 24,25,26].Despite this, the role of PYGO2 in MM remains unknown. Here, we evaluated the expression profile of PYGO2 in primary PCs of 1q21+ MM patients, and we characterized its functional roles in MM in vitro.MethodsPurification of primary CD138+ BM PCsWe analyzed a total of 72 purified BM CD138+ PCs from a cohort of MM patients. Of the 72 patients, 29 (11 SMM and 18 NDMM) are part of our previously published GSE227907 database [27]. The remaining patients, used as a validation set, are 24 NDMM and 19 RRMM, with a median age of 75 years and a range of 41–89 (Table S1). All participants in the study provided their written consent in accordance with the Declaration of Helsinki. The Ethics Committee of the Hospital of Parma approved all study protocols.Mononuclear cells (MNCs) were obtained from BM aspirates using Ficoll-Hypaque density sedimentation (Bichrome AG, Berlin, Germany). CD138+ MM cells were isolated by the immunomagnetic method using anti-CD138 monoclonal antibody-coated microbeads (MACS, Miltenyi Biotech, Bergisch-Gladbach, Germany). Only samples with >90% purity (flow cytometry analysis) were used.CoMMpass database analysesPublicly available data from the Multiple Myeloma Research Foundation CoMMpass Study (release IA22) were downloaded in March 2025 from the Research Gateway portal (https://research.themmrf.org/). Baseline BM aspirates with both whole-genome CNA segmentation files and RNA-seq gene counts files, as well as complete clinical annotations, were retained, yielding 847 unique patients. All genomic coordinates refer to GRCh38.Raw segmented log2 copy number ratios were corrected for purity using ABSOLUTE v1.0.6 [28] and for baseline region bias using BOBAFIT v1.0 [29]. Gene coordinates of genes relevant for MM biology and prognosis located on chromosome 1q (PYGO2, CKS1B, BCL9 and MCL1) were intersected with segments to derive the genes' CN deviations from baseline (CN_dev). Continuous CN_dev values were converted to discrete calls (deletion ≤ −0.2, neutral > −0.2 and < 0.2, gain ≥ 0.2, amp ≥ 1.2). Spearman correlations between CN_dev values of 1q genes and expression values were calculated.Raw RNA-seq gene counts provided by the CoMMpass consortium were normalised with the variance stabilising transformation in DESeq2 v1.42 [30] using the vst function. For each gene, quartiles were computed; the top quartile (Q1) defined high expression, while Q2–Q4 were combined as the reference group.Overall survival (OS) and progression-free survival (PFS) were measured from diagnosis to death or progression; patients without an event were censored at last follow-up. Kaplan–Meier curves were produced using survival v3.5.5 and survminer v0.4.9 R packages and compared with two-sided log-rank tests. Multivariable Cox proportional hazards models were used to compute Hazard Ratios (HR), reported with 95 percent confidence intervals.Binary matrices denoting common cytogenetic and mutational events were compiled: gain 1q (CKS1B), amp 1q, del 17p (TP53), del 13q (RB1), del 1p (CDKN2C), gain 8q (MYC), hyperdiploidy, t(4;14), t(11;14) and pathogenic TP53 single-nucleotide variants or indels. Frequencies in the PYGO2 Q1 group were compared with those in Q2–Q4 using Fisher's exact tests. Resulting p values were adjusted by the Benjamini–Hochberg method; a false discovery rate below 0.05 was considered significant. All analyses were executed in R v 4.3.2.Cancer dependency analysis (DEPMAP portal)Gene dependency of cell lines, gene expression, and copy number data for PYGO2 were obtained from the DepMap Portal Public 24Q4 release, downloaded in May 2025. The dependency data set used was CRISPR knockout gene effect scores computed with the Chronos model.MM cell lines were identified by using the DepMap cell line selector, by setting the Disease Subtype filter to “Plasma Cell Myeloma”. This selection yielded a total of 34 MM cell lines, while the remaining 1159 tumor cell lines served as the non-MM reference cohort.Spearman's rank correlation coefficients were calculated between PYGO2 Chronos scores, expression and copy number values in the MM cell lines. A one-tailed Wilcoxon rank-sum test was used to compare the PYGO2 Chronos scores between MM cell lines and other cancer types, as we hypothesized that MM cell lines have lower dependency scores compared to other cancer types. A two-tailed Wilcoxon rank-sum test was used to assess whether there is a significant difference in expression levels between MM cell lines and other cancer types. All analyses were conducted in R v4.3.2.PYGO2 Knock-down and PYGO2 OverexpressionPYGO2 knockdown in JJN3, OPM2 and U266 was obtained by infection with a pGFP-C-shLenti lentiviral vector anti-PYGO2 scramble (Scrb) and anti-PYGO2 (shPYGO2) (TL317702B, Origene). The sequence used for shPYGO2 was: 5’ CCTGCATACTCACATCTGACGGAGTTTGC 3’. PYGO2 overexpression (OE) in OCI-MY5 was obtained by infection with ORF lentiviral vector, LV[Exp]-CMV>hPYGO2[NM_138300.4]:IRES:mCherry:T2A:Puro (241205-1107qbzVectorBuilder). Recombinant lentivirus was produced by transient transfection of 293T cells following the standard protocol. Viral supernatant was collected after transfection, concentrated by centrifugation at 26,000 rpm, and resuspended in cell growth medium. JJN3, OPM2, U266 and OCI-MY5 cells were seeded at 1 × 106 cells/1 mL. All cell lines underwent the same puromycin selection (2 μg/ml) to minimize potential selection-related effects. The efficiency of the infection was assessed by flow cytometry using the GFP percentage for JJN3, OPM2 and U266, whereas mCherry was used for OCI-MY5. Stably transfected HMCLs were maintained in RPMI medium containing 20% FBS with 2 μg/ml puromycin until use.RNA-sequencing and bioinformatics analysisUniversal Plus mRNA-Seq with NuQuant kit (Tecan Genomics, Redwood City, CA) has been used for library preparation following the manufacturer’s instructions (library type: fr-secondstrand). RNA samples were quantified and quality tested by Agilent 2100 Bioanalyzer RNA assay or TapeStation RNA assay (Agilent Technologies, Santa Clara, CA). Final libraries were checked with Qubit 3.0 Fluorometer (Invitrogen, Carlsbad, CA) and Agilent Bioanalyzer DNA assay.Libraries were prepared for sequencing and sequenced on paired-end 150 bp mode on NovaSeq 6000 (Illumina, San Diego, CA).The number of reads (in millions) produced for each sample is listed in the table below.Raw sequencing data were processed following a standard bioinformatics pipeline. Base calling, demultiplexing, and adapter masking were performed using Illumina BCL Convert v3.9.31. Adapter sequences were masked during demultiplexing, with masked regions converted to “N” characters and the corresponding base quality scores set to 2 to facilitate subsequent trimming. Quality filtering and adapter trimming were carried out using ERNE. Reads were aligned to the human reference genome (hg38) using STAR3 (default parameters), a splice-aware aligner optimized for RNA-Seq data that identifies splice junctions and maps reads accordingly. Transcript quantification was performed with StringTie4 (default parameters), which reconstructs full-length transcripts, including multiple splicing isoforms, and provides gene- and transcript-level abundance estimates.Differential expression analysis was conducted using DESeq2, which fits a Generalized Linear Model (GLM) for each gene to compare expression levels across conditions. DESeq2 employs shrinkage estimation for dispersions and fold changes, enhancing the stability and interpretability of the results. This approach allows for a more quantitative assessment of differential expression, focusing on effect size rather than mere presence. Normalization was performed using the median-of-ratios method, and statistical significance was assessed using the Wald test.Regarding Gene Sets Enrichment Analysis (GSEA), we downloaded the Molecular Signatures Database (MsigDB). To ensure the credibility of the analysis results, we selected 10,000 permutations in the software. A significant Hallmark gene set was screened based on normalized enrichment score (NES), while screening was based on normalized p-value