Human Papillomavirus Deregulates the Response of a Cellular Network Comprising of Chemotactic and Proinflammatory Genes

Despite the presence of intracellular pathogen recognition receptors that allow infected cells to attract the immune system, undifferentiated keratinocytes (KCs) are the main targets for latent infection with high-risk human papilloma viruses (hrHPVs). HPV infections are transient but on average last for more than one year suggesting that HPV has developed means to evade host immunity. To understand how HPV persists, we studied the innate immune response of undifferentiated human KCs harboring episomal copies of HPV16 and 18 by genome-wide expression profiling. Our data showed that the expression of the different virus-sensing receptors was not affected by the presence of HPV. Poly(I:C) stimulation of the viral RNA receptors TLR3, PKR, MDA5 and RIG-I, the latter of which indirectly senses viral DNA through non-self RNA polymerase III transcripts, showed dampening in downstream signalling of these receptors by HPVs. Many of the genes downregulated in HPV-positive KCs involved components of the antigen presenting pathway, the inflammasome, the production of antivirals, pro-inflammatory and chemotactic cytokines, and components downstream of activated pathogen receptors. Notably, gene and/or protein interaction analysis revealed the downregulation of a network of genes that was strongly interconnected by IL-1β, a crucial cytokine to activate adaptive immunity. In summary, our comprehensive expression profiling approach revealed that HPV16 and 18 coordinate a broad deregulation of the keratinocyte's inflammatory response, and contributes to the understanding of virus persistence.


Introduction
Cervical cancer is the second most common cancer in women worldwide. More than 520,000 women are diagnosed with invasive cervical cancer each year [1]. Cervical and other anogenital carcinomas arise as result of an uncontrolled persistent infection with a high-risk type human papillomavirus (HPV), in particular types HPV16 and HPV18 [2,3]. A detectable cervicovaginal HPV infection in young women is close to 1-2 years [4] before it is cleared, suggesting that HPV can evade host immunity. Indeed, the infection cycle of HPV is one in which viral replication and release is not associated with overt inflammation [5,6] and HPV-specific adaptive immune responses are often weak or lacking in patients with progressive HPV infections [7][8][9][10].
Expression of specific viral oncoproteins, E6 and E7, is required for maintaining the malignant growth of cervical cancer cells [23]. To understand how HPV infection may alter KCs and evade PRR activation, direct protein interactions including the binding of the HPV E6 oncoprotein to IRF3 have been studied [24,25]. An OncoChip expression study showed that retrovirally expressed E6 and E7 efficiently downregulated type I IFN responses in keratinocytes, but surprisingly also upregulated the expression of pro-inflammatory cytokines [26]. Another early microarray study described downregulation of interferon-inducible genes in KCs containing episomal HPV type 31 [27]. These studies indicated that HPV-derived proteins could meddle with host immunity but the full spectrum of interference is within the limitations of these studies not visible.
We aimed at understanding the effects of high-risk HPVs on the immune response in KCs. First, we confirmed expression of the viral RNA receptors in undifferentiated and differentiated cells, while DNA sensor TLR9 was restricted to differentiated cells, and showed that HPV does not interfere with expression levels of the PRRs. Next, we focused our studies on undifferentiated KCs, since these are the target cells for latent infection with HPV. We generated expression profiles of several different control KCs and KCs harboring episomal copies of entire HPV16 or 18 genomes [28,29] on microarrays representing 24,500 well-annotated transcripts to study differences in the baseline gene expression by the presence of HPV. In addition, we studied differences in response to triggering the viral RNA PRRs with the synthetic dsRNA poly(I:C). Although HPV is a DNA virus, non-self dsDNA can serve as template for transcription into dsRNA by polymerase III and induce type I interferon and NF-Kappa-B through the RIG-I pathway [30][31][32]. Here, we show that HPVs were able to dampen a network of genes associated with activation of the adaptive immune response encoding antimicrobial molecules, chemotactic and pro-inflammatory cytokines, and proteins that are involved in antigen presentation, and that most of them are interconnected via IL1B.

Ethics statement
The use of discarded human foreskin, cervical and vaginal keratinocyte tissues to develop cell lines for these studies was approved by the Institutional Review Board at the Pennsylvania State University College of Medicine and by the Institutional Review Board at Pinnacle Health Hospitals. The Medical Ethical Committee of the Leiden University Medical Center approved the human tissue sections (healthy foreskin, healthy cervix, HPV16-or 18-positive cervical neoplasias) used for staining. All sections and cell lines were derived from discarded tissues and de-identified, therefore no informed consent was necessary.

Cell culture
Human epidermal KCs were isolated from foreskin, vagina, or cervix of unrelated donors [33] and established on a layer of lethally 137 Cs-irradiated mouse 3T3 fibroblasts. Passage 4-5 of primary KCs -devoid of contaminating cells -were grown in serum-free medium (Defined KSFM, Invitrogen, Breda, The Netherlands). Partial differentiation was induced by 1.8 mM Ca2+ for 24 hrs, terminal differentiation by placing KCs in single-cell suspension into serum-free medium containing 1.75% methylcellulose and 1.8 mM Ca2+ for 24 hrs [33]. KC cell lines maintaining episomal copies of HPV16 and HPV18 were created via an electroporation technique described previously [28,29] but without antibiotic selection. The cell lines were 100% HPVpositive. Southern analyses confirmed the recircularization and subsequent maintenance of episomal viral genomes at approximately 50-100 copies per cell (data not shown). The HPV-positive lines growed at similar rates with population doubling times of ,2 days) and, when placed in raft culture, all underwent the late stages of the virus life cycle, such as genome amplification, late gene expression, and virus production (data not shown). HPV-positive cells were grown in monolayer culture using E medium in the presence of mitomycin C-treated 3T3 fibroblasts [28,29] for passage 6-7, and adapted to serum-free medium for one passage before experimentation. All cells used were tested and found free of mycoplasm. Where indicated, cells were stimulated with poly(I:C) (25 mg/ml, InvivoGen, San Diego, USA). CCL5 and IL-1B concentrations in supernatants were determined using the Quantikine ELISA kits (R&D Systems, Minneapolis, USA).

Total RNA isolation and quantitative RT-PCR
Total RNA was isolated using TRIzol (Invitrogen, Breda, The Netherlands) followed by the RNeasy Mini Protocol (Qiagen, Venlo, The Netherlands). Total RNA (0.2 mg) was reverse transcribed using SuperScript III (Invitrogen) and oligo dT primers (Promega, Madison, USA). Triplicate PCR reactions were performed with 20 pmol of gene-specific primers and Taq DNA polymerase (Promega) using PCR conditions and primers as described previously for TLRs [34] and SPRR2A [35]. Predesigned primers and probe mixes for TLR3, CCL5, IL1B, RNASE7, NLRP2, and GAPDH were from Applied Biosystems (Foster City, USA). Threshold cycle numbers (Ct) were determined with 7900HT Fast Real-Time PCR System (Applied Biosystems) and the relative quantities of mRNA per sample were calculated using the DDCt method with GAPDH as the calibrator gene. The relative levels of mRNA were determined by setting the mRNA expression level of the lowest expressing control KCs to 1, unless otherwise indicated.

Analysis of differential gene expression
We fitted a linear model in limma v2.14.7 [39] with 'virus' (HPV-positive) and 'stimulation' (4 and 24 hrs) effects. We used a nested variable within 'virus' for the individual cell lines, where HVKp1 and HVK16 were the reference cells for the HPVnegative and HPV-positive groups, respectively. Multiple-testing corrected p-values [40] and log2 fold changes were extracted for different contrasts. For Table S1, the 4 and 24 hrs timepoints were combined into one F-test in limma. One-dimensional hierarchical clustering of log2 fold changes derived from limma was done in Spotfire DecisionSite 9.1 v19.1.977 using correlation as similarity measure and complete linkage.

Functional genomics analyses
Functional annotation of the groups of co-regulated genes identified by hierarchical clustering was performed using Anni 2.0 [41]. We used GenMAPP v2.1 [42] to overlay expression on the TLR signaling pathway, which was based on automatic extraction from KEGG [43] hsa04620 (7/17/09) with improved layout using PathVisio v1.0 beta software [44]. The edited pathway is available from GenMAPP and WikiPathways [45].
We used CORE_TF (www.lgtc.nl/CORE_TF) based on Trans-Fac 11.2 and Ensembl 49 [46] to identify over-represented transcription factor binding sites in promoters compared to a random set of 2966 promoters (1000 bp upstream+exon 1). Microarray probe EntrezGene IDs were converted to Ensembl Gene IDs using IDconverter [47], entries resulting in multiple or missing Ensembl Gene IDs were removed. The match cutoff was set to minimize the sum of false positives and false negatives; position weight matrices with a p-value for over-representation #0.01 and a frequency below 50% in the random set were selected.
The network was constructed using Ingenuity Pathways Analysis (IPA 7.6; IngenuityH Systems, Inc., www.ingenuity. com). The 663 HPV signature genes were filtered for the more extreme log fold changes to obtain a gene signature strongly affected by HPVs, and to get the number of genes below 500, which is the maximum limit of IPA for making a network. Genes not connected were deleted, the remaining HPV signature genes that were initially excluded as stated above were included to generate the final network consisting of 212 connected genes. All edges are supported by at least one reference from the literature, from a textbook, or from canonical information stored in the Ingenuity Pathways Knowledge Base.

Expression of viral pathogen recognition receptors in KCs
We determined the mRNA expression of Toll-like receptors and retinoic acid-inducible gene I (RIG-I)-like receptors in undifferentiated, partially and fully differentiated KCs. Expression of the small proline-rich protein 2A (SPRR2A) was used as a molecular marker of KC differentiation (Fig. S1A). Undifferentiated KCs were found to express TLR1, TLR2, TLR3, TLR5, TLR6, TLR10, RIG-I and MDA5 (Fig. 1A, 1B). Among the viral PRRs, TLR7, TLR8 and TLR9 were not detectable while TLR3, RIG-I and MDA5 were expressed. In parallel experiments, transcripts of TLR4 and TLR7-9 were readily detected in mRNA samples from Ramos B-cells and monocytes (Fig. S1B). The expression in KCs is largely in line with previous reports by others [13]. HPV-positive KCs showed essentially the same pattern of PRR expression (Fig. 1A, 1B). Additionally, real-time RT-PCR showed similar levels of TLR3 in HPV-negative and HPV-positive KCs (Fig. 1C). Upon differentiation KCs also expressed the DNA sensor TLR9, which was confirmed by immunohistochemistry in human foreskin and cervical epithelia (Fig. S2). TLR9 was also expressed in the differentiated layers of HPV-positive cervical epithelial neoplasias (Fig. S3). The absence of TLR4 expression in differentiated KCs, which was confirmed by expression microarray (see below), is consistent with work by others showing that TLR4 was only found in HaCat cells, but not in primary human KCs [16,48]. The pattern of TLR expression in differentiated HPV-positive KCs was similar to that in HPV-negative cells. Thus, HPVs did not affect mRNA expression of the tested PRRs.

HPV signature genes
We subsequently studied whether HPVs affected the signalling of PRRs using genome-wide expression profiling. Control KCs (n = 4) and KCs with episomal HPV16 or HPV18 genomes (n = 4) of foreskin, vaginal or cervical origin from eight different individuals were used to include biological variation. Since HPVs infect basal KCs, we focused on the viral PRRs expressed in undifferentiated cells, including TLR3, RIG-I and MDA5, which respond to the synthetic dsRNA poly(I:C) [13]. In agreement with the RT-PCR data, the presence of HPV did not change the expression of these PRRs (Table S1).
To obtain a robust signature of genes affected by HPVs, we selected differentially expressed genes between HPV-positive andnegative KCs at 0, 4 or 24 hrs of poly(I:C) stimulation with a false discovery rate (FDR) of 0.05 (1529 probes). Furthermore, we applied an absolute log2-fold change filter $1 to select genes that were at least two-fold up-or downregulated (663 probes representing 634 unique genes), designated ''HPV signature genes'' (union of genes in Venn diagram Fig. 2A, Table S2). The majority of HPV-specific differentially expressed genes were shared between all three (213) or two (150) conditions, with most overlap between 0 and 4 hrs. Notably, 219 genes were changed in the virus-positive group only after 24 hrs of poly(I:C) stimulation, showing that the effect of HPVs was more pronounced after poly(I:C) stimulation.

Poly(I:C) response in control KCs
We first focused on the effect of poly(I:C) stimulation in control KCs. While after 4 hrs (Fig. 2B left) we found 123 differentially expressed probes that were mainly upregulated, the response was more balanced and involved over 700 genes after 24 hrs of stimulation (Fig. 2B right). Many genes were upregulated, including pathogen-sensing receptors (RIG-I, MDA5, PKR), adaptor molecules (MYD88, TICAM1/TRIF, TICAM2/TRAM), and interferon regulatory factors (IRF1, IRF6, IRF7), see Table S1. These results are similar to a previous report showing that poly(I:C) stimulation induces antiviral and inflammatory responses in KCs [13]. Overlay of differential expression after 24 hrs of poly(I:C) stimulation on the TLR signaling pathway (KEGG hsa04620) showed upregulation of the Jak-STAT signaling pathway, triggered by temporary upregulation of IFNB1 after 4 hrs poly(I:C) through the TRAF3/TBK1 signal transduction route, resulting in upregulation of STAT1 and chemotactic cytokines CXCL10 and CXCL11. In addition, via TRAF6 the NFkappa-B signaling pathway was triggered, activating cytokines/ chemokines TNF, IL1B, IL6, IL8, CCL3, CCL4, and CCL5 (Fig.  S4). The cytoplasmic RNA sensing receptors MDA5 and RIG-I, which are not shown in the TLR signaling pathway, initiate signaling pathways that differ in their initial steps from TLR3 signaling, but converge in the activation of TBK1 and NFKB [13,49].

Deregulation of poly(I:C) response in HPV-positive KCs
The differentially expressed genes in the HPV-positive cells upon poly(I:C) stimulation largely overlapped with those in control KCs (Fig. 2B). Next, we studied the effect of the virus in the context of the TLR signaling pathway. Activation of the TLR signaling pathway in HPV-positive KCs upon 24 hrs of poly(I:C) stimulation was largely similar to the response in control cells (Fig.  S5). However, when directly comparing HPV-positive andnegative cells after 24 hrs of stimulation, relative downregulation of the adaptor TICAM1 and several cytokines (IL1B, IL6, CCL5/ RANTES) was evident. These results suggest that the dsRNA PRR signaling pathway is less activated in HPV-positive cells (Fig. S6).

Co-regulated genes downregulated by HPVs
We extended our analyses to the full set of HPV signature genes, and identified genes with similar expression patterns over the sample groups by unsupervised clustering (Fig. 2C, Table S2). The gene dendrogram was cut at six clusters to generate profiles of co-regulated genes (Fig. 2C, 2D). To identify transcription factors possibly involved in the coordinated expression changes, we analyzed the promoter sequences of the genes in each of these clusters for enrichment of predicted transcription factor binding sites [46].
The first three clusters contained genes that were downregulated in HPV-positive compared to HPV-negative cells. Binding sites for early growth response (EGR) family transcription factors, involved in differentiation and mitogenesis, were significantly enriched in these clusters (Table S3) , TICAM1 and RNASE7 show the HPV effect as well as the biological variation inherent to using KCs derived from different individuals and different tissues, combined with two different HPV types (Fig. 3A). Downregulation of CCL5 and TICAM1 was confirmed by qRT-PCR ( Fig. 3B and 3D), and ELISA showed lower CCL5 secretion in HPV-positive KCs upon poly(I:C) stimulation (Fig. 3C). For the small number of cluster 3 genes (15 probes), including pro-inflammatory cytokines (IL1B, IL1A, IL6), baseline expression (most likely activated by serum components) and upregulation at 4 and 24 hrs of poly(I:C)stimulation were suppressed in HPV-positive cells. These genes were already upregulated after 4 hrs of stimulation, and showed promoter enrichment of binding sites for Rel/NFKB family members and STAT5 (Table S3).
Interestingly, the majority of expression cluster 2 and 3 genes followed a similar pattern of suppressed poly(I:C) response, suggesting that many of these genes are downstream targets of PRR signaling. We focused on the antimicrobial molecule RNASE7, a member of the RNase A superfamily with broad-spectrum antimicrobial activity and ribonuclease activity [51,52], which was not known to be affected by viral infection. qRT-PCR confirmed RNASE7 upregulation upon poly(I:C) stimulation in control KCs, and suppression of poly(I:C)mediated upregulation in the presence of HPVs (Fig. 4A). Normal cervical epithelial cells expressed RNASE7 throughout the epithelia, and high expression was observed in the basal layer, the in vivo equivalent to undifferentiated KCs (Fig. 4B). In contrast, RNASE7 protein was not expressed in any of the layers of undifferentiated cells within a representative HPV-induced CIN3 lesion. These data suggest that by suppressing the gene activation of antimicrobial molecules such as RNASE7, HPVs evaded the innate antiviral responses of the host.

Co-regulated genes upregulated by HPVs
Clusters 4-6 contained genes that were specifically upregulated in the HPV-positive compared to HPV-negative cells. Cluster 4 genes (167 probes) included heat-shock response genes, cell cycle regulators and genes involved in replication initiation, transcription and splicing. These HPV-activated genes were downregulated upon poly(I:C) stimulation, but not to the same level as in control KCs. Binding sites for MEF2A, involved in the activation of stressinduced genes, and E2F, a family of transcription factors with a crucial role in the control of cell cycle that is indirectly activated by HPV E7, were enriched (Table S3). Cluster 5 (112 probes) contained cancer-related genes including tumor-promoting cytokines/chemokines and their receptors, e.g. CXCR7, of which the expression was higher in HPV-positive KCs irrespective of poly(I:C) stimulation. Many transcription factor binding sites were enriched, including motifs binding the oncoprotein MYC (Table  S3). Finally, the smallest cluster 6 (11 probes) included several antiviral response genes (TRIM5, ZC3HAV1, IFIT2, RARRES3, CXCL16) that were stronger upregulated in HPV-positive than in control KCs. Enriched binding sites included IFN-stimulated response element (ISRE), bound by transcription factor ISGF-3, and binding sites bound by interferon-response factors (IRFs).
In summary, the presence of episomal HPVs caused downregulation of genes involved in innate and adaptive immune responses as well as KC differentiation, while upregulated genes were involved in cell cycle, RNA and DNA metabolism. Overall, these data showed that HPVs induced coordinated changes in KC gene expression, detectable in unstimulated 'baseline' cells (mainly expression clusters 1, 5, majority of cluster 4) or after poly(I:C) stimulation (mainly expression clusters 2, 3, 6).

HPVs deregulate cellular networks
Understanding the network topology of gene and/or protein interactions may identify highly interconnected gene ''hubs'' targeted by HPVs. Therefore, we explored connections among the HPV signature genes based on literature and high-throughput database information collected in Ingenuity Pathways Analysis [53]. On the resulting network of 212 genes, we overlaid the expression log2-fold changes of HPV-positive versus control KCs after 24 hrs of poly(I:C) stimulation (Fig. 5). The center of the network was formed by the most interconnected gene IL1B, necessary for activation of the adaptive immune response [54], and IL6. IL1B and IL6 were downregulated, and connected to genes encoding cytokines and antigen presentation molecules that were also lower expressed in HPV-positive cells. We studied IL1B in more detail, since it represented a central target for HPVmediated suppression of both the innate and adaptive immune responses of KCs. RT-PCR data validated the microarray data showing that both the baseline and PRR-stimulated levels of IL1B were downregulated in HPV-positive KCs compared to control cells (Fig. 6A). Also, both the baseline and PRR-stimulated IL-1b secretion was lower in HPV-positive KCs (Fig. 6B). Secretion of IL-1b requires activity of both the TLR/NF-kappa-B and the inflammasome pathways [55]. The TLR/NF-kappa-B pathway activates pro-IL-1b expression, which is cleaved to active IL-1b by the inflammasome. In addition to the downregulation of pro-IL-1b, HPVs specifically downregulated the genes encoding inflammasome components NLRP2 in three of the four HPV-positive lines (Fig. 6C) and PYCARD/ASC, but not NALP3, possibly contributing to the observed lower level of IL-1b. The most interconnected upregulated gene of the network was CDKN2A, involved in cell cycle progression. Thus, by targeting highly interconnected genes, HPVs reprogrammed the gene network of KCs in favor of immune escape and cell proliferation of HPVpositive cells.

Discussion
We studied systematic differences in genome-wide expression profiles of control and HPV-positive undifferentiated (basal) KCs focusing on immune-related effects. The parallel analysis of several control and HPV16-and 18-positive KCs from several genital tissues ensured that the results can be generalized. The HPVpositive KCs expressed the full array of HPV genes and mimic latent HPV infection in vivo, which is also reflected by the fact that these cells display the entire differentiation-dependent HPV life cycle upon culture in organotypic raft cultures [28,29]. Our studies revealed that while KCs are well equipped to respond to viral pathogens, latent infection with HPV results in suppression downstream of the PRRs as reflected by lower expression levels of effector molecules involved in innate and adaptive immune response.
No difference was observed in expression levels of viral RNA PRRs TLR3, TLR9, RIG-I, MDA5 and PKR between control and HPV-positive KCs. We found that viral DNA PRR TLR9 was lacking in the basal layers in stratified squamous epithelia, but expressed in the suprabasal layers of the non-neoplastic epithelium. Previous studies suggested that E6/E7 expression affected  neither the expression nor the function of TLR9 [17], whereas others reported that E6/E7 expression resulted in loss of TLR9 expression [12]. Our data showed that forced differentiation of HPV-positive KCs resulted in the expression of TLR9, however, as HPVs inhibits differentiation this may appear as TLR9 loss similar to what was seen previously [12]. Thus, TLR9 is absent in the cells targeted by HPV, but other viral PRRs are expressed, including RIG-I that has been shown to indirectly function as a PRR for DNA viruses [30][31][32], suggesting that in essence undifferentiated KCs can sense HPV infection.
As there were no overt differences in the expression levels of PRRs, we focused on the interference of HPVs with the downstream pathogen-sensing machinery. First, our data showed that HPVs downregulated genes that have a direct antimicrobial function. Moreover, the presence of HPVs was associated with the downregulation of an array of pro-inflammatory and chemotactic cytokines, and antigen-processing and presenting molecules, and IL-1b and IL6 were the hubs in the center of this HPV signature gene network. Notably, the expression level of most of these genes was already lower at baseline. Poly(I:C), which triggers viral PRRs including TLR3 and importantly also RIG-I, increased their expression level in HPV-positive KCs albeit not to the same level as in control KCs. Previously it was shown that HPV31-positive KCs responded less well to interferon stimulation [27] and this fits with our own data showing that interferon-inducible genes (cluster 2) are downregulated. Apparently, this is not the only immune signaling pathway that is downregulated by HPV as our data reveal that also the TLR and the RIG-I-like receptor signaling pathways are suppressed in HPV-positive KCs. Notably, the failure of HPV31-positive KCs to respond to interferon was associated with downregulation of STAT1 (25). Specific downregulation of STAT1 was found only in our HPV16-positive KCs  (data not shown) suggesting that there may be a number of typespecific interactions with the host's immune system. Together these data suggest that HPVs dampen but do not block PRR signaling, and imply that the attraction of innate immune cells to the site of HPV infection, the subsequent initiation of adaptive immunity as well as the recognition of HPV-infected KCs is slowed down but not prevented. This clearly corresponds with the fact that it may take months or even a year to control HPV infections [4], and the increase in HPV-infected subjects capable of mounting an HPV-specific immune response in time [56]. Furthermore, it fits with the detection of HPV-specific memory responses after infection [9,57,58].
In particular, we found that HPVs downregulated toll-like receptor adaptor molecule 1 (TICAM1), a critical molecule in the TLR3 pathway that mediates NF-kappa-B and interferonregulatory factor (IRF) activation via downstream molecules TRAF3, TRAF6 and RIP1 [59]. Notably, the other poly(I:C) recognizing PRRs also malfunction in HPV-positive KCs suggesting that HPVs affect the TBK1 and NF-kappa-B signaling pathways downstream of the PRRs and implying that downregulation of TICAM1 is just part of the immune evasion strategy of HPVs. This is also illustrated by our finding that HPVs downregulated inflammasome components -needed to convert pro-IL-1b to the active form of IL-1b [60] -contributing to the lower secretion of IL-1b by HPV-positive cells. Of all candidate downstream targets IRF1 [25], IRF3 [24], the coactivator CPB [61], the IkB kinase complex [62], and the interferon-stimulated gene factor 3 (ISGF3) transcription complex [63] have been named as targets for either E6 and E7 proteins of HPV responsible for downregulating NF-kappa-B and TBK1 signaling. Others, however, have shown that E6 -instead of downregulating -may promote NF-kappa-B signaling [26,64]. Importantly, all of these studies relied on the overexpression of either one or both oncoproteins, which is more relevant for our understanding of HPV-transformed cells. The strength of our study lies in the use of KCs with episomal expression of the full array of HPV genes reflecting latent infection [28,29]. It would be of great importance to perform a genome-wide study of HPV-positive KCs during differentiation and interaction with (innate) immune cells thereby closely mimicking the situation in situ, but such an experiment would be technically challenging.
Non-cleared infection with high-risk HPVs leads to cervical and other anogenital carcinomas in which the virus genome integrates in the host genome [2,3]. The replication cycle of the virus is tightly coupled to the differentiation of basal KCs to stratified squamous epithelia and it is well known that HPVs inhibit KC differentiation [11]. In our expression data, this was reflected by concerted upregulation of cell cycle regulators and DNA/RNA synthesis, and downregulation of epidermis development and KC differentiation genes. CDKN2A, a critical cell cycle regulator upregulated by HPVs, was identified as one of the highlyconnected hub genes in the network of HPV signature genes. Similar results were described by Nees et al. using a cDNA oncochip [26].
We have shown that HPV16 and 18 dampen a cellular immune-related network in HPV-positive KCs, and affect a much broader spectrum of PRR responses than the previously described IRF route. Our study provides a framework for future exploration into the molecular mechanisms involved in HPV-downregulated immunity. The biological variation in gene expression between different donors might reflect genomic variation that could play a role the balance between clearance and persistence of HPV. Additionally, it would be of interest to study if other viruses capable of causing persistent infection or low-risk HPVs that cause benign genital warts use similar mechanisms to escape host's immune responses.  Figure S4 TLR signalling in KCs. Toll-like receptor signalling pathway (KEGG hsa4620) overlaid with differentially expressed genes between 24 hrs poly(I:C) stimulated and unstimulated uninfected keratinocyte cultures. Differentially expressed genes (FDR#0.05) were colored bright red (log2 fold change$1) or dim red (log2 fold change between 0 and 1) for upregulation upon poly(I:C) stimulation, or bright green (log2 fold change#21) or dim green (log2 fold change between 0 and 21) for downregulation. Grey boxes represent genes not fulfilling the above criteria, while white boxes are genes not represented by probes on the array. (PDF) Figure S5 TLR signalling in HPV-KCs. Toll-like receptor signalling pathway (KEGG hsa4620) overlaid with differentially expressed genes between 24 hrs poly(I:C) stimulated and unstimulated HPV-infected keratinocyte cultures. For explanation of colors, see Figure S4. (PDF) Figure S6 Differential TLR signalling between HPV-KCs and KCs. Toll-like receptor signalling pathway (KEGG hsa4620) overlaid with differentially expressed genes between HPV-infected and uninfected keratinocytes, both after 24 hrs poly(I:C) stimulation. Differentially expressed genes (FDR#0.05) were colored according to their log2 fold change (see legend Figure S4) for upregulation (red) or downregulation (green) in HPV-positive cells. (PDF)

Supporting Information
Table S1 Differential expression of pattern recognition receptors and signalling molecules in HPV-infected and uninfected keratinocytes. (PDF)