Introduction

Actinic Keratoses (AKs) are dysplasias of the epidermis that are associated with an increased risk for squamous cell carcinoma (SCC)1. It was estimated that in 2005 over 58 million people in the United States carred a diagnosis of AK2. AKs are associated with significant health care expenditure, with some studies estimating costs of over $1 billion per year (92% for physician office visits, 5% for prescription medications)3,4.

Untreated, AKs can either remain stable, regress, or progress to become squamous cell carcinoma (SCC). It is difficult to predict whether AK will progress to SCC, but there is a greater than six-fold increase in the risk of developing skin cancers with the presence of AKs5,6. It has been reported that SCCs originate in a lesion previously diagnosed as an AK with frequencies ranging between 65 and 97%4. In contrast, the estimated rate of progression from AK to SCC is less than 5%7. Additionally, AKs have been shown to be markers of field cancerization in excised basal cell carcinoma (BCC), SCC, and malignant melanoma (MM) specimens8.

The purpose of treating AKs is two-fold—to decrease risk of progression to SCC and to improve cosmetic appearance. Topical imiquimod is a field treatment for AK shown to lead to initial clearance rates of 85%. Sustained field clearance at 12 months was observed in 73% of patients treated with imiquimod9. Imiquimod’s antitumoral activity results from activation of the innate immune system and stimulation of antigen-presenting cells such as dendritic cells and macrophages. Imiquimod is an agonist for toll-like receptor 7, which is commonly involved in pathogen recognition and leads to downstream activation of nuclear factor kappa B and induction of pro-inflammatory Th1 cytokines including interferon alpha, tumor necrosis factor alpha, and various interleukins10,11,12,13,14.

Using the nanoString assay, we profiled the gene expression patterns in patients before and after treatment with imiquimod. Consistent with its function, we find that imiquimod diminishes gene expression associated with oncogenesis. We also find that imiquimod decreases expression of genes with immune suppressive function. Further, pre-treatment samples from complete responders have a distinct pattern of gene expression relative to partial and non-responders, suggesting that the pre-existing immune micro-environment in the AK may determine imiquimod response.

Methods

Tissue collection and initial histological examination were conducted at Icahn School of Medicine at Mount Sinai. The study was approved by the Insitutional Review Board (IRB) at Mount Sinai School of Medicine and was conducted in accordance with relevant guidelines and regulations (https://clinicaltrials.gov/ct2/show/NCT03914417). Informed consent was obtained from all participants and research was performed in accordance with the Declaration of Helsinki. Patients with at least four to eight visible AKs on the face and/or scalp were enrolled after informed consent was obtained. One AK from each patient was biopsied (3 mm punch) at day 0 and sent to pathology for diagnostic confirmation and future transcriptomic analyses. The remaining AKs were photographically documented, and one AK was listed as the target lesion. 21 patients with presumptive AK were enrolled in the study. After histological examination, 2 biopsies were found to contain seborrheic keratosis (SK), and not AK. These patients were therefore excluded from analysis, leaving 19 patients with analyzable pre-treatment biopsies.

Treatment with imiquimod 3.75% was initiated on a short cyclical treatment regimen. Imiquimod was applied daily, 2 weeks on/2 weeks off, for a maximum of 6 weeks (4 weeks of therapy with 2 weeks break) as described in a previous study15. The number of AKs was recorded at each follow-up visit. Biopsy of the target lesion was performed at week 14 if it was still present. If the target lesion was no longer visible, a biopsy was done at the site where it was previously located. Our study included 14 patients for whom matched pre and post biopsies were available, and 5 patients with unmatched biopsies (either pre biopsies or post biopsies, but not both).

In accordance with previously established methods, efficacy of treatment was measured by comparing the maximum lesion count during treatment (Lmax) to the number of AKs at week 14, two months after completion of treatment16. This allows for measurement of subclinical lesions that are unmasked after treatment with imiquimod. Number of AKs at Lmax and at week 14 was compared using a two-tailed paired t-test. Patients with AKs at the end of the trial were considered to be incomplete responders (IR) and patients without AKs at end of trial were considered to be complete responders (CR).

Pathology review, transcriptomic analyses, and statistical analysis were conducted at Columbia University Irving Medical Center with IRB approval (AAAO2758). All experiments were conducted in accordance with institutional guidelines and regulations. Formalin-fixed paraffin-embedded (FFPE) AK specimen blocks were measured and cut to provide a total of 250 mm2 of tissue. RNA was extracted using the miRNeasy FFPE kit (Qiagen) following the manufacturer’s protocol. Extracted RNA was quantitated by Agilent Bioanalyzer with RNA Nano chip assay.

RNA expression levels were profiled using a modified human PanCancer Immune Profiling panel (nanoString), which measures the expression of 770 genes plus 18 additional immune genes that were spiked into the panel, for a total of 788 genes17. RNA samples that passed quality and concentration standards were hybridized to target-specific probes and controls in a single tube for 20 h at 65 °C using 100–400 ng of RNA. Target-probe complexes were purified and immobilized on the nCounter prep station. Using the nCounter detection analyzer (nanoString), digital counts for each target RNA were acquired. Finally, nSolver software (nanoString) was used for normalization using housekeeping genes and for paired comparisons of gene expression. The NanoString assay includes RNA spike ins, labeled A-F in decreasing order of concentration, with positive spike in F (POS_F) in the raw data accounting for the lower limit of detection. Thus, transcripts with prohibitively low copy number are excluded (normalized to 0)18.

In order to determine whether relevant immune cell populations were in fact present in the AKs, we histologically examined the immune infiltrate present in five untreated AKs. 4 micron slides from pre-treatment biopsies were cut and stained with H&E, CD3 (T lymphocytes; clone LN10; Leica; 1:200 dilution), CD20 (B cells; clone L26; Thermo Fisher; 1:1000 dilution), or CD68 (macrophages; clone KP1; Biogenex; 1:100 dilution).

p values comparing the differential gene expression was calculated using the false discovery rate (FDR) correction. Gene expression in pre-treatment AKs was compared between two groups: (1) IR versus CR, (2) all patients with AEs compared to patients without AEs. Pathway analysis was performed as previously published19. Hierarchical clustering of genes in heat map was performed using a one minus Pearson correlation.

Given our intial set of genes, we applied pathFinder, a previously developed efficient graphical algorithm20, to extract a network from the background ConsensusPathDB (CPDB) signaling pathways. The core of pathFinder is the classical Depth First Search (DFS) algorithm, which expands on the initial input gene set by using genes located in the paths connecting input genes in the background network. Since the background network contains directed and undirected edges, in order to make all the edges in background network directed, each undirected edge was transformed into two directed edges with the same two end nodes but opposite directions. These two edges were not allowed to appear simultaneously in one path.

For every gene in the input list, the DFS explored all paths in the background network that started at that gene. The exploration of a path was stopped if it reached length k (here k = 2 was used) or arrived at a node with no valid child node(s). All nodes along the paths were included in the pathFinder output.

For Key Driver Analysis, we used the R package KDA21 (KDA R package version 0.1, available at http://research.mssm.edu/multiscalenetwork/Resources.html). The package first defines a background sub-network by looking for a neighborhood k-step away from each node in the target gene list in the network. Then, stemming from each node in this sub-network, it assesses the enrichment in its k-step (k varies from 1 to k = 3) downstream neighborhood for the target gene list. This algorithm takes the network structure and the location of input genes into consideration and allows us to discover key drivers significantly impacting genes in the input list.

Results

Imiquimod significantly decreases number of AKs

19 patients with confirmed histological diagnosis of AK were included in the final analysis (Table 1). The average age was 73 years and one patient was immunosuppressed, although not a transplant recipient. All patients had a history of prior AK and the treated lesions were located on the scalp (n = 11) and face (n = 8). AEs likely related to use of imiquimod occurred in eight patients, of which six had available pre-treatment biopsies. These AEs included but were not limited to pruritus, severe localized skin reaction, and facial swelling (complete list of AEs included in Table 1). Comparison of Lmax to number of AKs at end of study showed a significant decrease in AKs (t = 7.1, p < 0.0001, Fig. 1). Additionally, there was a significant decrease in the total number of AKs by the endpoint of the study at week 14 (t = 6.234, p < 0.0001, data not shown).

Table 1 Demographic information and treatment response in all patients included in the analysis. *Patients 15 and 17 initial biopsy showed SK and thus are excluded from this table.
Full size table
Figure 1
figure 1

Comparison of maximum number of AKs (Lmax), during treatment with imiquimod (black bars) with number of AKs at week 14, two months post-treatment with imiquimod (grey bars). Numbers above bars indicate number of AKs. Patients 15 and 17 had seborrheic keratoses on initial biopsy and are excluded from this figure. The samples with both pre- and post- treatment biopsies were: 1, 4, 5, 7, 8, 9, 10, 11, 12, 13, 14, 16, 17, 18, 19.

Full size image

Inflammatory genes are downregulated in post-treatment biopsies

A total of 32 RNA samples passed quality control checks and were analyzed using nSolver. NanoString results were filtered to search for genes with a change between pre and post-treatment samples (Fig. 2A). Levels of CDK1, CXCL13, IL1B, GADPH, TTK, ILF3, EWSR1, BIRC5, PLAUR, ISG20, and C1QBP were significantly higher in pre-treatment samples than in post-treatment samples (p < 0.05; Supplementary Table S1), with CXCL13 and IL1B showing a greater than 2 log2 fold decrease post-treatment. CXCL13 and IL1B are both implicated in both cancer metastasis and recruitment of immune cell populations to tumor cells, particularly B cells and macrophages22,23,24. The general function of these genes, as well as their role in oncogenesis, is described in Supplementary Table S2. Histological examination of five untreated AKs showed significant immune infiltrate of CD3, few CD20 positive cells, and moderate staining for CD68 in untreated AKs (Fig. 2B–E). Of the 16 post-treatment biopsies, 14 were characterized as scar tissue and 2 had residual AK. Thus, post-treatment AKs were not histologically examined. These results suggest that imiquimod diminishes expression of both immune modulatory genes and genes implicated in oncogenic processes.

Figure 2
figure 2

(A) Volcano plot of gene expression changes in post versus baseline of pre-treatment samples shows significant downregulation of CXCL13, IL1β, GAPDH, CDK1, TTK, ILF3, EWSR1, PLAUR (p < 0.05). (B) H&E 10 × magnification of pre-treatment AK shows dense immune infiltration (C) IHC staining of CD3 cells at 10 × magnification in pre-treatment AK shows significant staining for CD3 cells (D) IHC staining of CD20 cells at 10 × magnification in pre-treatment AK shows few CD20 positive cells (E) IHC staining of CD68 cells at 10 × magnification in pre-treatment AK shows moderate staining for CD68 cells.

Full size image

Complete responders have a distinct pattern of immune gene expression relative to incomplete responders at baseline

We next sought to identify differences in gene expression between IR and CR using the 788 gene panel. 103 genes were differentially expressed in pre-treatment AKs between IR versus CR, 95 of which were upregulated in CR (p values and genes shown in Supplementary Table S3). This was significantly more than was expected out of 788 total genes tested (p < 0.001, Fischer’s exact test). Further unsupervised clustering of these 103 genes distinguished CRs (patients 5, 1, 7, 11, and 9) from IRs (patients 19, 10, 14, 13, 12, 8, 4, 18, and 16), as shown using a heatmap in Fig. 3. Pathway analyses found that 15 pathways were significantly differentially expressed between CR and IR (Supplementary Table S4).

Figure 3
figure 3

Heatmap of 103 genes differentially expressed between complete responders and incomplete responders. Data shown is pre-treatment gene expression for genes that are differentially expressed between CR and IR. Blue indicates a lower expression, red indicates a higher expression. Heatmap was generated using Morpheus v.1 Software. URL: https://software.broadinstitute.org/morpheus/.

Full size image

Patients with adverse effects are more likely to respond favorably to treatment and exhibit a distinct pattern of gene expression

We next examined potential predicters of clinical outcome, including presence of AEs, history of SCC, history of smoking, immunosuppression, treatment area, and age. We found that presence of AEs significantly modulated treatment response, with AEs associated with a more favorable response to treatment (Fig. 4). The top ten most variable genes with respect to AEs are shown in Supplementary Table S5. We then compared patients with available pre-treatment biopsies who had AEs (n = 6) to those who didn’t have AEs (n = 10), finding 104 differentially expressed genes (Supplementary Table S6). 7 pathways were significantly different when comparing patients with AEs to those with no AEs (Supplementary Table S7).

Figure 4
figure 4

(A) Violin plot examining the effects of potential confounders on clinical outcomes. (B) Variance explained by gene for the first ten genes in the analysis.

Full size image

Network analysis shows densely clustered pathways involved in CR versus IR comparison

Finally, we generated a network from the differentially expressed genes found when comparing CR to IR, patients with AEs to patients without AEs, and CR with AEs to IR without AEs (see methods, Supplementary Figure S1). When comparing these networks, we found that the network utilizing the differentially expressed genes between CR and IR is significantly denser than the other two networks, suggesting that distinct network compartments are enriched by the two sets of differentially expressed genes. To confirm that these differentially expressed genes play a central role in highlighting different biological processes that regulates drug response and adverse effects, we performed pathway analysis. 350 pathways were found to be enriched by the network of CR vs IR, associated with different responses to treatment (Supplementary Table S8) and 4 pathways were found to be associated with the presence of AEs (Supplementary Table S9). 3 pathways overlapped between these two groups (type II interferon signaling, RANKL-RANK, and IL-4).

Discussion

Approximately 65% of SCCs arise from AKs7. Given that imiquimod is a first-line treatment for all AKs, it is important to better understand its effect on gene expression and the immune microenvironment. Further, there has been interest in determining whether gene expression levels can identify high-risk AKs, allowing for treatment of those lesions and not their indolent counterparts. Finally, it has been previously shown that the transcriptomic profile of AKs is more inflammatory and tumorigenic than normal skin25, and a better understanding of these differences may help in identifying high-risk AKs.

Our results further support the efficacy of imiquimod for treatment of AKs, given the significant reduction in number of AKs post-treatment relative to Lmax. This is in line with previous studies which have found that imiquimod 3.75% daily application using a short cyclical treatment regimen has similar efficacy as treatment with imiquimod 5% utilizing less frequent application over a longer treatment schedule15. We found that presence of AEs including pruritus, facial swelling, apthous ulcer, and flu-like symptoms was associated with improved treatment response.

Our NanoString findings highlight differential regulation of several genes associated with oncogenic pathways after treatment with imiquimod. Housekeeping genes are included in Supplementary Table S10. The genes that were significantly downregulated after treatment with imiquimod (CXCL13, IL1B, GAPDH, TTK, ILF3, EWSR1, BIRC5, PLAUR, ISG20, C1QBP, S100A and CDK1) are all associated with pro-oncogenic and immunoreregulatory pathways22,24,26,27,28,29,30,31,32,33,34,35,36. CXCL13, a chemokine ligand for CXCR5 that is associated with B-cell recruitment, may be involved in cancer progression and has been shown to be overexpressed in oral SCC37,38. IL1B, GAPDH, and CDK1 in particular have previously been shown to be overexpressed in various SCCs, including lung and head and neck SCC39,40. Interestingly, EWSR1 has been shown to inhibit the p53/p21 pathway involved in tumor suppression. Interruption of p53 is also associated with development of cutaneous SCC41. As such, the decreased expression of these genes post treatment with imiquimod is likely a reflection of the pre-SCC nature of the AK and of the disappearance of the lesion after treatment with imiquimod. Additionally, imiquimod eliminates the immunosuppressive tumor microenvironment.

Most importantly, we also find differences in pre-treatment gene expression levels between complete and incomplete responders to imiquimod. Of the 103 genes that were differentially expressed between the groups, the majority were upregulated in CR compared to IR (95 upregulated, 8 downregulated).

Further, pathway analysis showed that fifteen pathways were significantly different between CRs and IRs. Several pathways involved in interleukin signaling were significantly different between the CR and IR groups. These include the interleukin 1 receptor pathway, interleukin 6 signaling, interleukin 10 signaling, type 1 interferon and inflammatory cytokines production, and nuclear factor kappa B (NFKB)-related pathways. This suggests that inflammation and cytokine production may predict treatment response to imiquimod. The NFKB pathway has previously been shown to be a therapeutic target for head and neck SCC42.

Key driver analyses found three pathways that overlapped as significant key drivers for both response and presence of AEs. These included interleukin 4, receptor of the NFKB ligand (RANKL-RANK) signaling pathway, and type II interferon signaling pathways, all of which are immune-related pathways RANKL-RANK in particular is dysregulated in many types of cancer, including oral SCC, and may be a therapeutic target43. The exact role of these pathways on the evolution of AKs in response to imiquimod should be further elucidated using basic science techniques.

Additionally, our results show that presence of AEs may be associated with favorable response to treatment. This is important clinical information, as it highlights the importance of continuing treatment with imiquimod if mild, tolerable AEs occur.

Importantly, the post-treatment biopsies were performed two months after the final treatment with imiquimod. As a result, only 2 of the 16 post-treatment biopsies had residual AK, while the others were characterized as scar tissue. Thus, these gene expression changes are more likely to show the differences between AK and treated skin. Previous studies have examined gene expression changes during treatment with imiquimod and have found upregulation of inflammatory cytokines (IFN-alpha, IL10 receptor 1, TLR7) during treatment44. Additionally, overexpression of oncogenic genes and decreased expression of tumor suppressor genes has been observed in both sun-exposed, non-lesional skin and AKs. Treatment with imiquimod reverses this aberrant gene expression45. Another study performed biopsies of basal cell carcinomas treated with imiquimod as soon as the tumor showed signs of erosions and similarly found upregulation of genes involved in the inflammatory responses46. Our findings differ from these studies due to the later timing of the biopsy. Because the biopsy was performed after treatment and not during treatment, the inflammatory state had likely already resolved, allowing us to examine changes in gene expression between AK and treated skin.

Potential limitations of this research include the small cohort size and the lack of information about progression to SCC. Additionally, the majority of post-treatment biopsies were scar tissue. Analysis of gene expression in treatment-resistant AK may provide more information about resistance to imiquimod. Moreover, as discussed above, a biopsy during treatment with imiquimod may provide additional information regarding mechanisms of clearance.