At the outset of this work, no human m²G RNA methyltransferases had been identified. Therefore, to expand the repertoire of human RNA methyltransferases and to explore the roles of the enzymes and the relevant modification in human cells, we sought to identify human m²G RNA methyltransferases.

To this end, we performed homology comparisons to known m²G RNA methyltransferases of other species, and we identified human THUMPD2, THUMPD3 and TRMT11 as potential m²G RNA methyltransferases (Fig 1A and 1B). This highlighted TRMT11, THUMPD2 and THUMPD3 as prime candidates for being human m²G RNA methyltransferase. Indeed, during the course of our study, all three human proteins were confirmed as m²G RNA methyltransferases [19–21]. Of these, we also identified that THUMPD3 is more highly expressed (more than 2-fold) in lung adenocarcinoma (A549 and H1975) cells than in normal IMR-90 lung fibroblasts (Fig 1C). Given this elevated expression, we decided to further explore THUMPD3 in lung tumourigenesis.

In initial analyses, siRNA and shRNA approaches were employed to deplete THUMPD3 from human lung adenocarcinoma cells (A549). The use of two distinct methods to target THUMPD3 minimised the risk of following off-target effects. As previous studies had implicated THUMPD3 in the methylation of tRNA species [19, 20], we first confirmed that the depletion of THUMPD3 from A549 cells reduced m²G in tRNAs. We therefore performed RNA MS analysis on tRNA isolated from control and THUMPD3 depleted A549 cells. As expected, the results confirmed a significant decrease of m²G in the small (≤200 nucleotides) RNA fraction following either siRNA or shRNA mediated THUMPD3 depletion (Fig 1D and 1E). Collectively, these findings indicate that THUMPD3 expression is elevated in lung cancer cells where it functions as an m²G tRNA methyltransferase.

To determine whether THUMPD3 is required for lung cancer cell growth, we adopted an siRNA approach to deplete the protein from two NSCLC cell lines, A549 and H1975. Live-cell imagining confirmed that depletion of the enzyme from either lung adenocarcinoma cell line hindered their proliferation (Fig 2A and 2B). In addition, THUMPD3 depletion from A549 cells via an independent shRNA approach also had a negative effect on cell proliferation (S1A–S1C Fig). Crucially, the growth defect was completely rescued by exogenous expression of non-targetable THUMPD3 (Fig 2C). Importantly though, depletion of THUMPD3 from normal, non-transformed lung fibroblasts did not significantly influence their proliferation rate (S1D Fig), suggesting a cancer cell specific dependence of THUMPD3.

To further investigate the role of THUMPD3 in lung cancer, the colony formation capacity of A549 and H1975 cells was assessed. THUMPD3 depletion from either lung cancer cell line impaired their colony formation (Fig 2D and S1E Fig). In addition, we discovered that overexpression of THUMPD3 not only enhanced lung cancer cell proliferation, but it also stimulated the proliferation rate of normal lung fibroblasts (S1F and S1G Fig). This highlights a potential oncogenic role for THUMPD3 in human lung cells.

We next investigated whether depletion of THUMPD3 impacts the metastatic potential of lung adenocarcinoma cells, a crucial hallmark of cancer [22]. For these and the following experiments we predominantly used an siRNA approach to streamline experiment numbers, but importantly all of our previous analyses indicated no obvious differences between siRNA or shRNA mediated THUMPD3 depletion. The effects of THUMPD3 depletion on cell migration were firstly evaluated using wound healing (scratch) assays. Wound healing was monitored for 24 hours, using light microscopy, after a scratch was introduced. Notably, decreased wound healing capacity was evident in THUMPD3 depleted cells compared to control cells, as early as 8 hours after scratch initiation (Fig 2E and 2F). This effect was recapitulated in cells where shRNA targeting of THUMPD3 was used (S1H Fig). To further investigate the effect of THUMPD3 depletion upon cell migration, a standard Transwell assay was conducted. Cells with siRNA mediated depletion of THUMPD3 displayed reduced migration compared to cells transfected with control siRNA (Fig 2G). Overall, these data support a role for THUMPD3 in lung cancer cell migration.

To further investigate whether the cell proliferation defect following THUMPD3 depletion is associated with apoptosis, we examined apoptotic cell death mediated by caspases. Caspases execute apoptosis by cleaving several essential proteins crucial for cellular function and survival, with PARP-1 being one of the most well-known caspase substrates. Caspase-mediated cleavage of PARP-1 is widely recognized as a hallmark of apoptosis [23]. To assess apoptosis following THUMPD3 depletion, we utilised an antibody that specifically detects the 24 kDa fragment of cleaved PARP-1. This indicated elevated levels of cleaved PARP-1 72 hours post-knockdown of THUMPD3, suggesting that THUMPD3 depletion induces apoptosis in (S1I Fig).

The above experiments demonstrate that depletion of THUMPD3 hampers the proliferation and migration of lung adenocarcinoma cells. However, the underlining mechanisms of these changes remain unknown. To gain deeper insights into pathways responsible for these observed phenotypes, we conducted a whole-transcriptome analysis on A549 cells using RNA-sequencing (RNA-seq) technology. THUMPD3 was depleted using two different approaches: siRNA and shRNA, with cells harvested 5- and 6-days post treatment initiation, respectively.

PCA analysis of RNA-seq data indicated strong sample and data quality, as evidenced by the clear distinction between untreated and treated samples, while untreated samples themselves showed no significant variation (S2A Fig). Additionally, replicates of the same experimental condition clustered closely together, further confirming the consistency and reliability of the data.

Differential expression analysis identified a total of 284 transcripts that exhibited significant upregulation in both targeting approaches, while 189 gene transcripts were significantly downregulated in both approaches (Fig 3A, S2E–S2G Fig and S1 Table). Reassuringly, THUMPD3 was one of the most downregulated transcripts, confirming efficient mRNA targeting and depletion.

We then conducted Reactome enrichment analysis to gain insights into biological processes influenced by THUMPD3 depletion (Fig 3B, S2B and S2C Fig). One of the most upregulated pathways was associated with extracellular matrix (ECM) organisation. Amongst downregulated pathways, key processes such as cell cycle regulation were identified. Overall, these findings establish an important role for THUMPD3 in various biological pathways likely critical for the efficient propagation and survival of lung cancer cells.

Real-time quantitative polymerase chain reaction (qPCR) validation of four transcripts encoding ECM proteins (Fibrillin-1, Fibronectin-1, Tenascin C and Laminin B3) confirmed they were significantly upregulated upon depletion of THUMPD3 (Fig 3C). Surprisingly though, when we assessed the levels of the 3 ECM proteins encoded by the upregulated transcripts, we found that they were significantly reduced following THUMPD3 depletion (Fig 3D and S2D Fig). This suggests that despite elevated ECM transcript levels in THUMPD3 depleted lung cancer cells, either their translation or the stability of the relevant encoded proteins was impaired. Additionally, it is well established that cells tend to downregulate their ECM proteins upon loosing cancerous properties [24], as they do with THUMPD3. An intriguing question remains as to whether the observed reduction in ECM protein levels is a direct or indirect consequence of THUMPD3 depletion. Ultimately though, the changes we see in ECM proteins following THUMPD3 depletion likely underpin, at least in part, the impaired cellular migration phenotype in THUMPD3 depleted lung cancer cells.

Given the recently identified link between THUMPD2 and RNA splicing [20, 21], we explored whether THUMPD3 was also involved in this process. While RNA-seq can be used to evaluate overall steady-state transcript levels, the resulting data can also be interrogated for changes in the relative abundance of gene isoforms, with alternative splicing (AS) being one main cause. Using stringent criteria, we evaluated our dataset for potential AS events induced by THUMPD3 depletion. In total, we identified 137 alternatively spliced events (S2 Table). In terms of AS event classification, intron retention emerged as the most prevalent (40%), followed by alternative exon usage (25%; Fig 4A). We generated volcano plots for the primary AS event types: alternative exon usage, intron retention, alternative 3’ and 5’ splice sites (Fig 4B and S3A–S3C Fig). Despite intron retention being the predominant alternative splicing event type detected, a review of the literature revealed no clear connection between the identified targets and THUMPD3. Additionally, we found no obvious links between these transcripts and the phenotypic changes observed in lung cancer cells following THUMPD3 depletion. In contrast, however, transcripts exhibiting alternative exon usage events were enriched in mRNAs associated with ECM and cell adhesion molecules (CAMs; S3D Fig, Table 1). In addition, substantial number of identified targets had some links to neurodevelopment (Table 2).

Interestingly, one of the most statistically significant alternative exon usage events upon THUMPD3 depletion occurred within FN1 mRNA (encoding Fibronectin-1), one of the upregulated ECM-related transcripts (Fig 4A and S1 Table). THUMPD3 depletion induced a highly selective decrease in usage of exon 25 within FN1 (Fig 4C). Strikingly, exon 25 encodes a unique domain known as extra domain B (Fig 4D), which is associated with cancer-related functions [25, 26]. Importantly, qPCR analysis confirmed that THUMPD3 depletion results in exon 25 loss (Fig 4E).

Cell culture was carried out under sterile conditions in a standard laminar flow hood. Human cell lines were maintained in filter top flasks or culture dishes and cultured in respective growth media. 293T, A549 and IMR-90 cell lines were grown in DMEM supplemented with 1 x penicillin-streptomycin (pen-strep) antibiotics and 10% [v/v] foetal bovine serum (FBS). H1975 cells were grown in RPMI media supplemented with 1x pen-strep antibiotics and 10% [v/v] FBS. Cells were grown at 37°C, 5% [v/v] CO₂ and were maintained at ≈ 80% viability. Cells were passaged every 3–5 days and seeded at 1 x 10⁵ live cells per ml. Cells were passaged using conventional cell culture techniques. The maximum number of passages was 20 before a new vial of cells was revived. Tests for mycoplasma contamination were carried out every month or when a new vial of cells was revived.

ON-TARGETplus (Dharmacon) or FlexiTube (Qiagen) siRNAs were used for the knockdown experiments according to the manufacturer’s instructions. 0.5–3 x 10⁵ cells per well were used for transfections in 6-well plates. siRNA sequences are listed in S3 Table.

1 x 10⁵ cells per well were plated into 6-well plates. Cells were allowed to settle and attach overnight. Next morning, transfection was carried out using Lipofectamine^™ 3000 Transfection Reagent (Thermo Fisher, L3000001) according to the manufacturer’s instructions.

The ‘all-in-one’ pLKO-Tet-On system was used to generate inducible shRNA mediated THUMPD3 knockdowns [57, 58]. 293T cells were seeded onto L-polylysine (1 μg/ml) coated 10 cm dishes 24 hours prior to transfection in antibiotic-free DMEM media with 10% [v/v] FBS. Cells were 80% confluent on the day of transfection. To produce lentiviral particles, cells were transfected using 42 μl FuGENE^® 4K Transfection Reagent (Promega, E5911) with 6 μg of l pLKO-Tet-On constructs together with 3.5 μg psPAX2 (Addgene, 12260) and 4 μg pCMV-VSV-G (Addgene, 8454). 16 hours later, the media was exchanged for fresh media. 48 hours post-infection, virus particles were harvested and sterile filtered using 0.45 μm syringe filters (Millipore). Aliquots of lentiviral supernatants were stored at -70°C. Target cells (A549) were then transduced with the virus with polybrene (8 μg/ml) to increase the efficiency of transduction. 48 hours after infection transduced cells were selected for by treatment with puromycin (1 μg/ml). To induce expression of shRNAs, 10–100 ng/ml doxycycline was added to the media and the cells were incubated for 4–8 days to allow adequate expression of the shRNA. Target sequences of shRNA are listed in S3 Table.

0.5–2 x 10⁵ cells per well (A549, H1975 or IMR-90) were seeded in 6-well dishes and reverse transfected with 2.5 nM of control (siscr) and THUMPD3 (siD3 #3, 9) siRNAs. Next day the media was changed, and plates were placed into an Incucyte^® machine for live-cell imaging. The confluence of cells at time point zero was maintained at ∼ 20% to allow a proper curve to be formed. The algorithm is based on measuring cell proliferation using live-cell time-lapse imaging without labels using Classic Confluence Analysis.

Stable cell lines were developed as described by Ebrahimi et al., 2015 [59]. Essentially, the pCMV6-Entry (pEV) and pCMV6-THUMPD3 (pD3.3) (S3 Table) were linearised using BglII enzyme (NEB, R0144S) according to the manufacturer’s instructions. Following electrophoresis of the resulting digestion product in a 1% [w/v] agarose gel, the linear plasmid was purified using standard QIAquick Gel Extraction Kit (Qiagen, 28706X4). H1975 cells were transfected with 1 μg of the purified linearised vector as described above. Next day culture medium was refreshed. 48 hours post transfection, the cells were treated with 400 μg/ml of G418 antibiotic to select for transfected cells. Once control cells without plasmid died, the concentration of G418 was reduced to 175 μg/ml. Successful overexpression of THUMPD3 in the relevant stable cell lines was confirmed by Western blotting with anti-THUMPD3 antibody. Cell lines with equivalent expression of endogenous and exogenous THUMPD3 were selected for further analyses, such as rescue experiments.

H1975 cells stably harbouring control (pEV) or THUMPD3 (pD3.3) expression vectors were used. 3 x 10⁴ cells per well were seeded into a 12-well plate and reverse transfected with 2.5 nM of control (siscr) and THUMPD3 (siD3 #9) siRNAs. 24 hours later the media was exchanged, and the plate was placed into the Incucyte^® instrument to follow cell growth.

0.5 ml/well for 12-well or 1 ml/well for 6-well dish of Crystal violet (CV) solution (0.05% [w/v] Crystal violet, 20% [v/v] ethanol) were added to fix and stain cells. Cells were incubated in the solution for 5 minutes on a shaker (slow shake). Plates were washed by submerging in a 1 litre beaker with H₂O, allowed to dry overnight and imaged on an Epson Perfection V800 Photo scanner.

A549 (1000 cells/well) and H1975 (2000 cells/well) were seeded into 6-well plates in a total volume of 2 ml per well. Colonies were allowed to form over 10–14 days. Colonies were then fixed and stained with CV and imaged using an Epson Perfection V800 Photo scanner.

A549 or H1975 cells were reverse transfected with 2.5 nM of control (siscr) and THUMPD3 (siD3 #3, 9) siRNAs. 72 hour later 2 x 10⁴ cells were seeded into ibidi chambers (Culture-Inserts, 2 wells/chambers; ibiTreat, 80206) and allowed to form a monolayer. 24 hours later the Culture-Inserts were removed (’scratch’ was introduced). The chamber was filled with 1 ml of medium. To ensure that images were taken at the same field, a line was drawn perpendicular to the scratch at the bottom of the imaged field. Each image was then taken at the same spot. Pictures were taken at 0–24 hours post scratch using a Leica EC3 digital camera. Images were analysed using Image J software. The area of the wound was measured using the Image J plugin—MRI wound healing [60].

Transwell^® inserts (8 μm; Corning^®, 354480/1) were used to assess the migration capacity of cells upon THUMPD3 depletion. 2 x 10⁵ cells (A549 or H1975) were reverse transfected with 2.5 nM of control (siscr) and THUMPD3 (siD3 #3, 9) siRNAs. After 72 hours, cells were serum starved for 24 hours. The following day, 800 μl DMEM 10% [v/v] FBS were added to the bottom of the empty well (12-well plate). Transwell^® chambers were then added to the wells with sterile tweezers, ensuring no bubbles were introduced. The cells were washed, detached with 1x trypsin-EDTA or cell dissociation buffer enzyme-free PBS-based (Gibco, 13151–014) and neutralised in DMEM 2% [v/v] FBS. The cells were counted and prepared at a concentration of 1 x 10⁵ cells per ml in serum-free media. 300 ml of cell suspension (3 x 10⁴ cells) was added to each Transwell^® chamber. The plates were then incubated overnight at 37°C. 24 hours post-seeding, the Transwells^® were removed and washed by submerging in 1x PBS. The migrated cells were then fixed and stained using CV, and imaged on an Epson Perfection V800 Photo scanner.

Protein samples were harvested in RIPA buffer (25 mM Tris-HCl pH 7.6, 150 mM NaCl, 1% [v/v] NP-40, 1% [w/v] sodium deoxycholate, 0.1% [w/v] SDS) supplemented with protease and phosphatase inhibitor tablets, followed by sonication in Bioruptor Pico, Diagenode. Lysate was cleared by centrifugation at 12,000 × g and quantified using a DC Protein Assay Kit (BIO-RAD). Samples were run on SurePAGE^™ precast polyacrylamide gels and transferred to nitrocellulose membrane. Membranes were blocked in 5% [w/v] BSA (or powdered milk) in 1x TBS, 0.1% [v/v] Tween20 (TBS-T) and then incubated at 4°C with primary antibodies. Bound antibody was detected with ECL solution following the manufacturer’s instructions (Promega, W1015). The membrane was imaged via a Chemidoc^™ imaging system (BioRad). Antibodies and reagents listed in S3 Table.

RNeasy Mini Kit (QIAGEN, 74104) was used for RNA extraction. An additional step of DNase digestion was incorporated after the first wash following the RNase-Free DNase Set protocol (QIAGEN, 79254). To prepare cDNA, RNA was reverse transcribed using the SuperScript^™ III Reverse Transcriptase (Invitrogen^™, 8080044). qPCR analysis was performed using a StepOnePlus^™ Real-Time PCR System (Applied Biosystems^™, Thermo Fisher, 4376357) with Power SYBR Green Master Mix Kit (Applied Biosystems^™, Thermo Fisher, A25741). Primers listed in S3 Table.

RNA Clean & Concentrator kits (ZYMO RESEARCH, R1013 or R1017) were used to separate small (≤ 200nt) RNA fraction, following the manufacturer’s instructions. RNA concentration was determined using a Qubit^™ RNA HS Assay Kit (Thermo Fisher, Q32855).

Nucleosides were prepared from enzyme-processed RNA by enzymatic digestion, using a cocktail of Benzonase (Merck), Phosphodiesterase 1 (Merck), and Antarctic Phosphatase (New England Biolabs) as described previously [61]. The reactions were filtered using an Amicon 30kDa MWCO spin-column (Merck) to remove protein and the filtrate was mixed with a 2 x loading buffer containing 0.1% formic acid and an internal standard (13C-labeled uridine generated from 645672–1MG Merck KGaA, previously treated with Antarctic Phosphatase). The samples were loaded onto an ACQUITY UPLC HSS T3 Column, 100 Å, 1.8 μm, 1 mm X 100 mm (Waters Corp., Milford, MA, USA) and resolved using a gradient of 2%–10% acetonitrile in 0.1% formic acid over 10 min. MS analysis was performed in positive ion mode on an Orbitrap QExactive HF (Thermo Fisher, Waltham, MA, USA) mass spectrometer. Standard dilutions of all experimental nucleosides were prepared and analysed in parallel. There were three technical replicates of each sample and the analytical processing was performed using XCalibur Software (Thermo Fisher).

Whole transcriptome analysis was performed using RNA-sequencing (RNA-seq) technology. Total RNA, depleted of rRNA, was used as an input material for RNA-seq. rRNA was removed using Ribo-Zero^® rRNA Removal Kit (Illumina) following the reference guide. RNA Clean & Concentrator-5 kit (ZYMO RESEARCH, R1013) was used to clean up the rRNA depleted RNA sample. RNA was eluted in 15 μl nuclease-free water. The efficiency of rRNA removal was assessed using an RNA Screen Tape (Agilent, 5067–5576) analysed on a 4200 TapeStation System (Agilent, G2991BA). 13 μl of rRNA depleted RNA sample were used for the preparation of the RNA-seq library. The NEXTFLEX^® Rapid Directional RNA- Seq Kit 2.0 protocol was followed to perform all the remaining steps for library construction steps.

Reads were aligned to the human reference genome (GRCh38.p14) using STAR (v2.7.10a), with –quantMode GeneCounts to quantify gene counts with the annotation source NCBI RefSeq reference transcriptome GTF file [62]. Processing steps including gene filtering, normalisation and subsequent differential gene expression analysis were performed using the DESeq2 (v1.38.3) package in R (v4.2.0) [63].

In the DESeq2 workflow, two experimental variables were considered for the experimental design. The first variable was ‘Type,’ representing the type of gene knockdown (shD3 #1, shD3 #4, siD3 #3). The second variable was ‘Condition,’ representing untreated versus treated samples. Three biological replicates were used for each condition. Differential analyses were performed on all types of knockdowns using the formula: ∼ Condition + Type + Condition:Type. Additionally, analyses were conducted on each knockdown type individually using the formula: ∼ Condition, where treated samples were compared to untreated samples (reference level). Adjusted p-values were calculated using the Benjamini-Hochberg method as implemented in DESeq2 [64]. Fold changes were shrunk using the ashr and apeglm methods as implemented in DESeq2 [65, 66].

Functional enrichment analyses, encompassing Gene Ontology (GO), Reactome and KEGG pathway analyses, were conducted using the clusterProfiler package (v4.7.1.3) in R [67]. The list of significantly differentially expressed genes from DESeq2 was used for these analyses.

Alternative splicing (AS) events between untreated and treated shD3 #1 samples, representing conditions without and with THUMPD3 depletion, respectively, were identified from aligned reads using MAJIQ (Modelling Alternative Junction Inclusion Quantification) (v2.4) software [69, 69] and its associated visualisation package VOILA (v2.4). Default parameters were used, including –threshold 0.2 and –changing-pvalue-threshold 0.05. The AS information of Local Splicing Variations (LSVs) was further parsed and classified into subtypes of AS events (such as exon skipping, alternative intron, 5’ alternative splice sites, 3’ alternative splice sites etc.) using the MAJIQ Modulizer program with stringent parameters: 1) maximum dPSI value ≥ 0.2, 2) p value for determining whether a LSV is changing ≤ 0.05, 3) probability that a LSV is changing ≥ 0.95.