Non-invasive diagnosis of non-melanoma skin cancer subtypes predicted by microRNAs recovered from adhesive tape discs

We investigated whether a tape-stripping approach could be used for analysis of microRNA (miR) expression as a non-invasive method for prediction of non-melanoma skin cancer subtypes. From direct miR sequencing of RNA recovered from adhesive tape discs applied to skin lesions prior to biopsy, we determined rank lists of miRs to predict superficial or invasive BCC or SCC. Expression of these miRs was then assessed by digital PCR of disc-derived miR from independent validation cohorts of BCC (n=24 Sup, n=38 Inv) and SCC (n=45 SCC-in-situ, n=14 Inv). Expression ratios or products of >2 miRs were superior to individual miRs in discriminating Sup from Inv BCC and SCC subtypes. Moreover, several miR combinations distinguished Sup BCC and Inv SCC from a cohort of benign keratoses that can clinically resemble these cancers. Overall design: Cohort identification and RNA preparation for sequencing Prior to biopsy of lesions, an adhesive disc (D-Squame D101, Clinical and Derm) was firmly applied to the lesion multiple times with rotation to ensure the entire adhesive surface had made physical contact. Discs were immediately stored at -20 oC and then at -80 oC until processing. Whole RNA was isolated from the discs using a miRNeasy Micro Kit (217084, Qiagen) and the RNA eluate was stored at -80 oC until sequencing. Analysis of miR sequencing The human GRCh38 genome and gene annotation files were downloaded from Ensembl release 106 and a reference database was created using STAR version 2.7.9a (Dobin et al., 2013). A custom Perl script (smallRNA_pe_umi_extractor.pl from https://github.com/HuntsmanCancerInstitute/UMIScripts) was used to trim adapters and extract the Unique Molecular Index (UMI) from the second read of a paired-end Qiagen small RNA library with 76x51 base pair reads. The trimmed reads were aligned to the reference database using STAR options optimized for shorter miRNA reads (>=16 bases matched to the genome, <= 5% mismatches over mapped length, splicing switched off). A second Perl script (bam_umi_dedup.pl) removed PCR duplicates from the aligned BAM file and the mapped reads were assigned to annotated genes in Ensembl and to miRNAs in miRbase release 22.1 (Kozomara et al., 2019) using featureCounts version 1.6.3 (Liao et al., 2014). Differentially expressed miRNA in miRbase were identified using a 5% false discovery rate with DESeq2 version 1.36.0 (Love et al., 2014). Principal components analysis was run using the regularized log values from the top 500 variable miRNAs and ellipses were added representing the 90% confidence interval from a multivariate t-distribution. Statistical Method for Feature Selection Features (i.e. miRs) were identified as described previously (Fastner et al., 2024). Briefly, binomial lasso regression, as implemented in the R (version 4.3.2) package “glmnet” was used to select a group of features that separated subtypes (e.g. superficial from invasive BCC). Ten-fold cross validation was used to select the penalty parameter for the lasso regression. Low expression miRs with an adjusted mean count of less than 10 were removed from consideration. qRT-PCR assays RNA isolated from discs was converted to cDNA using a TaqMan™ Advanced miRNA cDNA Synthesis Kit (Thermo Fisher Scientific, A28007) and stored at -80 oC until used. Samples were read using QuantStudio™ Absolute Q™ Digital PCR System (Thermo Fisher Scientific) and analyzed using QuantStudio™ Absolute Q™ Digital PCR Software (version 6.3, Applied Biosystems) according to the manufacturer's recommendations. Samples were diluted 1:100 in 0.1X TE buffer (Ambion, AM9849) and digital PCR was performed using a MAP16 Digital PCR kit (Ambion, A52688) according to the manufacturer's instructions. Briefly, each reaction consisted of 2 ?L Absolute Q DNA Digital PCR 5X mix, 0.5 ?L 20X TaqMan Advanced miRNA probe (Applied Biosystems, A25576) and 7.5 ?L cDNA, with 15ul of Isolation Buffer added over each sample prior to sealing the wells with the provided gasket strips. Data was analyzed using Quant Studio Absolute Q digital PCR software, version 6. Threshold bars were set using the manual option and results reported as number of copies per ?L. Normalized expression values were obtained copy number for an individual miR by that for the invariant miR-16-5p. Combined expression scores were obtained from the products and/or quotients (based on whether miRs were “up” or “down,” respectively) of normalized expression scores of multiple miRs. Pearson correlations of normalized digital PCR values with normalized regularized log sequencing values was determined using R version 4.3.2. ROC curves and other analyses were generated using Prism software.