Mechanistic signatures of HPV insertions in cervical carcinomas

To identify new personal biomarkers for the improved diagnosis, prognosis and biological follow-up of human papillomavirus (HPV)-associated carcinomas, we developed a generic and comprehensive Capture-HPV method followed by Next Generation Sequencing (NGS). Starting from biopsies or circulating DNA samples, this Capture-NGS approach rapidly identifies the HPV genotype, HPV status (integrated, episomal or absence), the viral-host DNA junctions and the associated genome rearrangements. This analysis of 72 cervical carcinomas identified five HPV signatures. The first two signatures contain two hybrid chromosomal–HPV junctions whose orientations are co-linear (2J-COL) or non-linear (2J-NL), revealing two modes of viral integration associated with chromosomal deletion or amplification events, respectively. The third and fourth signatures exhibit 3–12 hybrid junctions, either clustered in one locus (MJ-CL) or scattered at distinct loci (MJ-SC) while the fifth signature consists of episomal HPV genomes (EPI). Cross analyses between the HPV signatures and the clinical and virological data reveal unexpected biased representation with respect to the HPV genotype, patient age and disease outcome, suggesting functional relevance(s) of this new classification. Overall, our findings establish a facile and comprehensive rational approach for the molecular detection of any HPV-associated carcinoma and definitive personalised sequence information to develop sensitive and specific biomarkers for each patient.


INTRODUCTION
Each year, 4500,000 individuals are diagnosed worldwide with human papillomavirus (HPV) associated cancers. 1 HPV is a hallmark of most cervical cancers but is also associated with anal and head and neck carcinomas. Upon persistent, long-term infections of the cervical mucosa, oncogenic HPVs can cause high-grade cervical intraepithelial neoplasias (CINs), that without ablative treatment, may develop into invasive carcinomas. 2,3 Better understanding of HPV biology has permitted major advances in the development of prophylactic vaccines and in the detection of HPV-associated diseases. However, cervical carcinomas remain frequent worldwide and their prognosis is poor, particularly in advanced stages and in cases of relapse. 4,5 Notably, the efficiency of surgical removal of tumour relapses is related to the tumour size. 6 It is thus of major importance to acquire molecular biomarkers that will facilitate the early diagnosis of cervical cancer, the detection of precocious tumour relapses and offer optimal biological follow-ups of the disease, upon and after treatment. Therefore, efforts to facilitate the diagnosis of HPV-associated cancers, elucidate the mechanisms driving HPV oncogenicity and develop cost-efficient and precise monitoring tools to improve clinical oncology, are needed.
At the molecular level, HPV contains a 7.8-kbp double-stranded DNA circular genome (Figure 1a). One major difficulty in HPV detection and characterisation of tumours lies in the highly polymorphic nature (over 200 genotypes and variants) of this virus and in the sporadic nature of viral integration in the host cell. A subset of HPV genotypes, able to infect the mucosa of the ano-genital and oropharyngeal tract, have oncogenic properties and are classified as low, medium and high risk, 7 among which the HPV16 and HPV18 genotypes are the most prevalent. [8][9][10][11] Viral oncogenicity is largely related to their encoded E6 and E7 proteins, which have the ability to interact and inhibit the host tumour suppressors, p53 and pRB, as well as many other cellular protein targets. 12 In CIN cases, the HPV DNA generally replicates extrachromosomally in the cell nucleus, whereas in most invasive carcinomas, a portion of the viral DNA is stably integrated into a unique site or perhaps within a few sites in the tumour cell genome. [13][14][15][16][17][18][19][20] Upon integration, there is often a disruption of the E1/E2 genes, which causes constitutive expression of the E6/E7 oncoproteins. Viral integration may also lead to structural and/or modified expression of cellular genes nearby, as well as genomic amplifications, encompassing the HPV insertion site that may have a role in tumour progression. [21][22][23][24] Along this line, the viral-cellular host junctions represent highly specific molecular markers of the individual patient tumours. 25 However, PCR-based approaches used for junction analysis 26 are difficult to implement in clinical oncology, notably due to a high rate of inconclusive and false-negative results. 23 Most recently, a few genome-wide studies based on Next-Generation Sequencing (NGS) techniques such as whole genome, exome and targeted capture as well as RNA sequencing 17,20,27 have been reported. However, these 1 approaches remain limited and/or are difficult to implement in clinical practice, owing to technical complexities, biases towards certain HPV types, thresholds of bioinformatic filtering and/or cost efficiencies.
Here we describe a sensitive and exhaustive double-Capture HPV method, followed by NGS, designed to determine any HPV genotype and to detect cellular integration sites in HPV-associated carcinomas. This approach has been validated from a panel of 72 cervical patient cases, using tumour biopsy and circulating DNA from blood samples. Surprisingly, our patient-to-patient analysis of the junction reads identified five cellular HPV status whose signatures reveal new insights into the mechanisms of HPV integration and concomitant formation of large deletions or duplications at the insertion sites. Furthermore, these distinct HPV integration signatures prove to be relevant, both biologically and clinically.

Capture-NGS of HPV DNA
The main steps of our Capture-HPV-NGS method are illustrated ( Figure 1 and Materials and Methods). Briefly, genomic DNA from frozen biopsies and blood samples were extracted and the HPV-containing regions enriched using a double-capture method. Both rounds of capture consist of hybridisation with an exhaustive set of 22,000 single-stranded biotinylated HPV probes (Supplementary Table S1) in order to enrich for the small targeted viral sequences, which represent less than 8 kb of DNA out of 46 Gb of human DNA (an~10 − 6 ratio). The probes were designed to detect the entire length of the HPV genomes from 4200 HPV genotypes or variants, allowing an exhaustive and unbiased method to capture any chromosomally integrated, or episomally non-integrated HPV DNA sequence. The HPV-enriched DNA library was then sequenced to obtain single-end long reads (~500 nt average length) using Roche 454 Life Sciences (Roche Company, Branford, CT, USA) Junior or paired-end reads (2 × 150 nucleotides (nts)) using Illumina MiSeq (Illumina Inc., San Diego, CA, USA; Materials and Methods). These read lengths allowed the characterisation of the HPV type and of the insertion junctions. In 2 out of the 50 HPV-integrated samples (#83 and 202), one of the HPV-chromosomal junction was unambiguously mapped but the second junction read could not be definitively oriented as it overlapped a repetitive region of the human genome. To be noticed, as a proof of concept workflow, we successfully performed barcoded multiplexing of up to 7 and 12 patient  Figure 1. Capture-NGS method to identify HPV genotypes and integration sites. (a) The HPV dsDNA circular genome contains early (E) and late (L) replicating genes, necessary for viral replication and capsid formation. Integration and disruption/loss of the E2 suppressor gene sequences may lead to overexpression of the E6 and E7 viral oncogenic sequences, as well as disruption of host genomic DNA sequences. (b) Capture-NGS of HPV-integrated viral-cellular junctions. Genomic DNA is mechanically and homogeneously nebulised, followed by dsDNA library preparation and annealing of adaptor barcodes. HPV-biotinylated single-stranded probes (~22,000) then bind and capture HPV-containing library sequences that are enriched on a streptavidin column. These enriched DNA products are sequenced by NGS whose read alignments permit identification of the genotype, as well as the HPV-cellular hybrid junctions and localisation of the HPV insert in the tumour genome. and HPV68 in two cases each (2.8%) and one case each of HPV6, HPV42, HPV51 and HPV52 (1.4%) ( Table 1). These figures confirm the ability of our approach to characterize the two most prevalent worldwide HPV16 and HPV18 viral genotypes, as well as the next less common or rare types which vary by geographic region, such as HPV31 and HPV33 in Europe, HPV31 and HPV45 in Americas or HPV52 in Eastern Asia, HPV52 in Eastern Asia or HPV45 and HPV73 in Oceania. 28,29 Of note, two cases of SCC associated with the low risk HPV6 and HPV42 genotypes were found. Importantly, even with the increase of sensitivity of the HPV capture, only one viral genotype was detected in each carcinoma case, demonstrating the clonal HPV status in each tumour. In two cases (No. 45, 46; Table 1), no HPV DNA sequences were detected upon mapping against the entire set of HPV variants using a permissive alignment threshold of 20 nts. Consistently, no specific PCR products were detected using a pair of HPV generic primers. 10 These tumours are definitive examples of HPV-negative cervical carcinomas.
Heterogeneity of HPV inserts. Looking at the 50 integrated cases (2J and MJ categories, Figure 2a; Table 1), we noted that the viral reads did not map uniformly along the entire HPV sequence but instead, several gaps identified as junction reads were observed, reflecting the presence of a deletion in the viral genome. Altogether, we found that a portion of the HPV genome, ranging from 25.2 to 96.8%, is inserted into the human genome (Supplementary Table S2). The breakpoints occur most often Abbreviations: EPI, episomal HPV genomes; HPV, human papillomavirus; NA, not available; ND, not determined; NEG, negative for HPV; rr, repetitive region; 2J-COL, two hybrid chromosomal-HPV junctions whose orientations are co-linear; 2J-NL, two hybrid chromosomal-HPV junctions whose orientations are non-linear. CGH array: flat/normal, deletion or amplification (FDA) structural profiles based on comparative genomic hybridisation (CGH) microarray analyses. Bold entries represent directly disrupted genes.

HPV integration signatures in cervical carcinomas A Holmes et al
within the early E1-E5 and late L1-L2 genes. All viral inserts contain the HPV origin of replication (ORI), located in the LCR region, suggesting that the process of integration may have positively selected for the retention of these elements. Most viral inserts contain disrupted or deleted E2-coding regions while retaining the E6 and E7 genes, a feature consistent with enhanced expression of the E6 and E7 oncogenes via the inactivation of the E2 repressor upon integration. However, unexpectedly, in one HPV51-associated case (No.6), the two breakpoints were located within the E6 and E7 oncogenes and no evidence of associated episomal HPV genomes was observed. Also unexpectedly, in two other cases (No. 12, 18), the E1 region was disrupted but the E2 gene remained intact, indicating that E2 disruption is not always associated with HPV integration. In three additional cases (No. 22,28,32), disruption of the E6/E7 oncogene was also found in a sub-set of the reads detected. These data outline the unappreciated and large heterogeneity of HPV sequence insertion patterns amongst different tumours, which validates the overall benefit of using our comprehensive and unbiased set of Capture-HPV probes.
Architectures of HPV insertions. To better understand the nature of the HPV integration events, we examined the viral-chromosome junction sequences and the integrated HPV content. The analysis of our NGS reads was performed on a patient-to-patient basis that provided a puzzling heterogeneity of HPV signatures and crucial mechanistic insights. This is reflected by the variable number of hybrid junctions captured per sample and their chromosomal distributions, being clustered at the same chromosomal locus or scattered at different loci. In terms of the number of breakpoints, we identified two hybrid junctions (2J) per tumour in 30/70 cases, indicating a simple integration pattern of a truncated viral genome into a single chromosomal locus (Table 1 and  Supplementary Table S2). Most interestingly, close examination of these 30 cases led us to focus on the relative orientation of the viral-chromosomal junction reads. As illustrated in Figure 2a, two classes of carcinomas, harbouring two hybrid junctions could be identified: In 11 cases (No. 1-8, 17, 126, 208), the upstream and downstream junction reads were in the expected orientations for a simple co-linear HPV integration event, i.e.,; in the orientation chromosome-HPV at the upstream breakpoint and HPVchromosome at the downstream breakpoint ( Two junctions, co-linear (2J-COL) Hereafter, these distinct viral insertion architectures, containing two junctions are referred to as 2J-colinear (2J-COL) and 2J nonlinear (2J-NL), respectively. In the remaining cases (No. 83, 202) two junctions were also identified but one of them being in a repeated region of the genome, their relative orientation could not be determined. In a third category (12 cases; No. 21-29, 115, 116, 144; Table 1, Figure 2c), the architecture of insertion exhibited 3-9 distinct hybrid viral-cellular junctions, each validated by at least five independent and non-overlapping reads, that were clustered in the same chromosomal region. This multiple junction and clustered pattern was abbreviated as MJ-CL. A fourth pattern ( Figure 2d) was characterised by the presence of several hybrid junctions (4)(5)(6)(7)(8)(9)(10)(11)(12) in which the chromosomal sequences aligned ectopically to unlinked and distant regions on the same chromosome or were scattered among 2 or more chromosomes. These cases are referred to as the multiple junction scattered (MJ-SC) category. On the basis of the analysis of the individual HPV polymorphisms, a single viral genotype was detected in all junction reads of the same tumour, again underlining the clonal and personalised HPV signature in each tumour. The final category (20 cases) includes high-coverage NGS reads containing full-length HPV genotypes and no hybrid reads represented by more than one NGS read. Therefore, this category, named EPI, corresponds to cases with only high copy episomal HPV genomes, without clonal integration (Table 1). To be noted, in one of these samples (#132), we detected 16 non-duplicate overlapping HPV-chromosomal reads (0.3%) suggesting that one copy of HPV was already integrated in a subset of the tumour cells. Our capacity to detect such heterogeneity is likely owing to the high sensitivity of the capture method that allows to reach very high sequencing coverage (up to several thousand folds  We then analysed homologies between the viral and cellular sequences at the breakpoint fusions. We detected no significant continuous homology lengths but rather microhomologies at and surrounding these breakpoint sites, an important feature, also described by others, 20 that should be considered when deciphering possible mechanism(s) of HPV viral integration (Supplementary Figure S2).
The distinct architectures of HPV insertions are also detected in circulating DNA. The above case studies utilised DNA from biopsy samples, allowing the identification of pertinent personalised genomic biomarkers. However, a tumour biopsy is an invasive procedure as compared with taking a blood sample. For the same reason, the analysis of circulating tumour DNA (ct-DNA) samples would also be preferable for the optimal follow-up of patients during and after treatment. Therefore, we investigated whether our Capture-HPV method could also detect the HPV integration signature biomarkers from blood samples. We analysed the ctDNA (Materials and Methods) of serum specimens taken at the time of diagnosis in five patients in which the 2J-COL (cases 5,7), 2J-NL (cases 18,20) and MJ-CL (case 29) signatures had been found, respectively. Capture-HPV from serum/plasma specimens was conducted as for biopsies (Materials and Methods) and yielded numerous HPV and viral-human hybrid NGS reads. Most importantly, the same viral-junction patterns were identified in the blood and tumour samples from each patient (Supplementary Figure S3). Furthermore, as for the biopsies, the coordinates of the paired junction reads delineate the associated chromosomal rearrangement at the HPV insertion site. Thus, the Capture-HPV approach is a powerful and robust method to detect HPV and characterise its genotype, physical status and integration loci, as well as the coordinates of the associated structural alterations. It should be emphasised, that the detection of ct-DNA by Capture-NGS does not require prior analysis of the matched tumour specimen to identify the desired tumour-specific molecular biomarkers.
Viral load and integration patterns. The viral load was measured with two different methods. Classically by quantitative PCR of the E7 gene copy number (Supplementary Table S2)  The highest number of reads (mean of 181,742) was observed in the EPI cases, a feature indicating a high copy episomal HPV content, consistent with the high viral load also detected by quantitative PCR (Supplementary Table S2, Supplementary Figure S4B).   Table 1, Supplementary Table  S2). Remarkably, they are often closer to each other in the 2J-COL (⩽5 kb in 6/11 cases) and more distant in the 2J-NL (⩽5 kb in 0/17 cases). These observations raise the important question of what happened to the recipient chromosomal region and whether this could provide clues on the mechanism of integration. To characterise the eventual associated rearrangements, we used a genome wide CGH-array approach that allows high-resolution (5 kb) determination of a copy-number change (Materials and Methods) implicating the region bracketed by the HPV insertion breakpoints identified by capture, and also detects chromosome gains/losses throughout the entire genome. We first examined the 2J-COL and 2J-NL categories. Concerning the 2J-COL events, in 5/11 cases, the distance between the HPV breakpoints were ⩽ 5 kb. In the other cases (No. 1-3, 5, 7) the chromosomal copy-number decreased about 1.5-fold between the two HPVcellular junctions (Figure 3a), indicating that the intervening chromosomal sequences between the breakpoints have been deleted, and that the insertion of HPV occurred on one of the chromosome pairs. In the first case (No. 2) we observed a 6.2-kb deletion directly within the OFD1 gene (Xp22.2). In the second case (No. 3), a 196.5-kb deletion is located near the ZBTB18 and C1ORF100 genes (1q44). Finally in the third case (No. 5) we found a 351.8-kb deletion leading to the partial deletion of both the     (Figure 3d). In addition, HPV DNA was often found to co-localise at the transition boundaries between lengthy chromosomal gains and losses (No. 33, chromosome band 20q11.21, see Figure 3d). The loss of heterozygosity observed in the vicinity of HPV insertion sites points to a mechanism of integration involving a break-induced replication type event. In sum, strikingly, no deletions were detected at HPV insertion sites in any category other than 2J-COL suggesting that, all of the other events result from the same 'ends-in' mechanism of chromosome sequence assimilation (see below). We also assessed by FrAGL analysis (Methods), the frequency of amplicon gain or loss and copy-number variations genome wide, among all cases. No additional striking genomic differences were observed among the signature groups (Supplementary Figure S7). Overall, the cases within the 2J-COL category share, as a common structural feature, the deletion of the chromosomal region within the interval of the two breakpoints. In contrast, the 2J-NL and MJ categories only show local amplifications. Amplification events at the HPV insertion locus have been previously described, 17,23,27 but to our knowledge, this is the first comprehensive analysis integrating the relative architecture of the two junctions in individual tumours and the detection of a deletion pattern arising from HPV integration. In contrast, to the usually large amplification events that are likely to affect several adjacent genes, the small size of the associated deletions are more likely to yield a single-gene dysfunction, more amenable for the design of precision therapeutics.

Annotations of the HPV insertions
As previously observed, the viral-host junction reads in HPV-associated carcinomas map at multiple locations in the human genome 15,17,20,23,30 Herein, we identified 175 robust viral-chromosomal sequence fusions (Supplementary Table S2 and Supplementary Figure S10). Among the 50 insertion cases, viral sequences were found within a coding sequence in 29 of these (58%) or within o500 kb from a known gene in 18 cases (36%; Table 1). In the remaining cases, HPV was inserted in integenic and repeated regions of the genome (Table 1). Recurrent integration events were found at the KLF5/KLF12 locus (five cases), and near MYC (three cases). Among the genes not previously reported in the literature as directly targeted by HPV integration, we found RB1, AKT3, SST and ID1, implied in the regulation of cell proliferation, LPP in cell adhesion/motility, transcription factors such as AFF3, BCL6, CCAT1 and CCAT2, and oncogenes such as RAB11A and RAB22A (Supplementary Figure S11). Among the potentially deregulated genes located near HPV DNA, we noted MAST4 and MAP2 implicated in the regulation of microtubules, MMP12 and COL4A4 coding for matrix components and PF4V1 in angiogenesis. Altogether, our data identify new cellular targets for viral integration in HPV-associated carcinomas and confirms that integration sites, though widely distributed, are not random, being most likely mechanistically favoured and/or functionally selected when carcinogenesis ensues.
Clinical annotations related to HPV integration signature In order to assess the relevance of our classification of cervical carcinomas, we further compared HPV signatures with biological and clinical parameters. No significant correlation was observed between the integration pattern, tumour size (17-80 mm), histology and clinical staging (I-IV). In contrast, a strikingly unbalanced distribution of HPV genotypes was observed according to the integration signatures. Unexpectedly, a majority (11/18) of the HPV18 associated cases were found in the 2J-NL group whereas only three of the most prevalent HPV16 genotypes (a total of 49% in the present series) exhibited this signature (Figure 4). In this category, we also detected one HPV45 case (that has 81.7% sequence identity with HPV18, Supplementary Figure S12) and two HPV73 cases. Also striking, we found heterogenous genotypes in the EPI category but no HPV18 nor HPV45.
Another aspect that we analysed concerns the age of the patients. The clinical data from our case panel gave a median patient age of 48 years. The median age was 48 in patients with 2J-COL tumours, 45 in 2J-NL patients, 45 in MJ-CL, 52 in MJ-SC and 56 in women with EPI carcinomas (Supplementary Figure S14). Altogether, an 8-11 year difference in median age was observed between the 2J and EPI classes (Pairwise comparisons using Mann-Whitney Wilcoxin rank sum test, 2J-COL versus EPI, P = 0.264; 2J-NL versus EPI, P = 0.099; MJ-CL versus EPI, P = 0.185; and MJ-SC versus EPI, 0.959) and between all the integrated versus EPI, the P value is 0.12. Finally, follow-up data recorded at least 6 months after diagnosis were available for 58 cervical cancer patients. Half were either deceased (14/29), diagnosed with pelvic relapse (10/29) or metastasis (3/29). By comparing the HPV signatures with the disease outcomes, we noted that a poor outcome was more frequently recorded in tumours from the 2J category (17/24) than in the EPI (9/17) and MJ (5/17) categories (Chi-square test: P value integrated versus EPI = 0.12; Supplementary Table S2). Further studies including a larger number of cases analysed in the present way are required to assay the robustness of this intriguing distribution for its potential prognostic value. The ability to determine the HPV mode of integration upon diagnosis of the pathology, will be of significant value for the follow up of the patient.

DISCUSSION
By studying tumour DNA and ct-DNA from cervical carcinoma patients, we validated a cutting-edge NGS method able to detect and characterise all tumour-associated HPV sequences. Our approach provides new insights into the mechanisms of HPV integration as well as its diverse structural consequences on the tumour cell genome at the insertion locus. Most importantly, this work has led to the classification of HPV-associated carcinomas into five categories that are potentially relevant both biologically and clinically. As such, HPV-associated tumours represent a highly favourable model for the development of genomic medicine.
Compared with PCR-based methods for HPV typing, the comprehensive Capture-HPV approach is unbiased and able to provide single-nucleotide level information, as deep genome-wide based-sequencing approaches, provide. 17,20,27 Further, our approach is rather cost efficient (few hundred dollars per sample), rapid (few days without the need to reiterate PCR attempts),  technically relatively simple for a molecular biologist, and is therefore able to be implemented today in clinical settings, with the benefits of providing in one experiment: (i) the capture of any HPV genotype regardless of its prevalence in the population, (ii) its full sequence identity by the de novo assembly of the NGS reads, (iii) the viral load assessment as measured by the read coverage normalised against house keeping gene(s) included in the capture probes, (iv) the physical status, integrated or episomal of the viral sequences, (vi) the pattern of integration into the diverse chromosomal insertion sites, (vii) the sequence of the hybrid viral-chromosomal junctions that provide exquisitely specific biomarkers, as well as (viii) the nature and coordinates of the associated chromosomal deletions or amplifications. This method also provides insights into the tumour biology and potential therapeutic opportunities to treat HPV-associated cancers, by characterizing the status of the HPV E6/E7 oncogenes and their E1/E2 repressor as well as all other possible viral rearrangements and the identification of the gene(s) that are directly or indirectly altered via HPV insertion.
In a pioneer study, the group of E. Schwarz developed a valuable approach combining ligation-mediated-PCR with NGS to specifically detect HPV16-associated sequences. 17 Integration of HPV16 DNA was observed in 78.8% of their cases, which correlates well with the 71.4% cases of integration found in our series, collected in two French hospitals. In comparison with previous approaches, one of the assets of Capt-HPV is the incorporation of a large set of probes representing 4200 HPV genotypes and variants that allowed us to detect and identify all HPV sequences, including those of various worldwide geographic prevalences such as HPV16, HPV18 and HPV45 in America, HPV45 and HPV33 in Africa or HPV52 in Eastern Asia, 28,29 as well as their status of integration. Furthermore, our case-by-case analysis of the HPV sequences led us to identify 5 HPV signatures distinguished by the relative orientation of the hybrid viral-cellular junctions (2J-COL versus 2J-NL), the proximity of HPV breakpoints in the 2J-COL category and their larger spacing in 2J-NL, the multiplicity of junction reads in certain tumours (MJ-CL, MJ-SC) and finally the presence of the viral genome in the episomal state only (EPI).  Figure 5. Ends-out and Ends-in HPV integration mechanisms. In the left panel (ends-out repair), the two linearised microhomologous 3' ends of the HPV fragment align co-linearly with the target chromosome and perform strand invasion by ends-out replacement of the chromosome sequence with HPV, resulting in a stable deletion of the intervening host sequences. In the right panel (ends-in repair), the two microhomologous 3′ ends align to the target chromosome sequences, only when HPV is in a gapped episomal form. The 3' ends of the gapped episome can invade and initiate synthesis from the target sequence, creating a chimeric viral chromosome. Repair synthesis occurs by microhomology mediated break-induced replication pathway (MM-BIR). The resulting viral-chromosome chimeric episome can then integrate by homologous recombination into the same copied target sites, resulting in chromosome duplications and amplifications, a transient and recurrent event, predicted to be unstable.
Remarkably, the distinction between the 2J-COL and 2J-NL integration patterns was independently validated by two characteristic genomic alterations found at their respective loci, namely deletions (2J-COL) or amplifications (2J-NL). These correlated features can be simultaneously explained by two modes of integration. In the first mode, the two ends of the linearised intact gapped viral DNA align to the host target, such that by 'ends-out' replacement, the virus is inserted creating a deletion of the intervening chromosome segment, explaining the formation of the 2J-COL architecture ( Figure 5, left panel). This is a model of DNA double-strand break repair proposed to explain plasmid/fragment integration in yeast. 31 Not surprisingly, this mechanism of integration that does not require extensive strand synthesis, favours closely spaced interactions of the HPV ends-out DNA ends, with nearby chromosomal regions. In the second mode of integration, the two ends of the linearised intact or gapped viral DNA, will invade the host target in an 'ends-in' mode, involving template directed DNA synthesis ( Figure 5, right panel), explaining the formation of the 2J-NL architecture and the opportunity for each HPV end to independently interact with distant regions of the chromosomal regions. This leads to the copy (via gap repair synthesis) of the target sequences at variable lengths, resulting in capture of the target sequence incorporated into the viral DNA episome, creating a chimeric viral-chromosome intermediate. Next, as part of an amplification 'looping' model 27 the initiating intermediate will be able to form DNA concatemers and/or integrate into the chromosomal insertion site by homologous recombination, thus explaining concomitant duplication(s)/ amplification(s) at the target sequence. In our seventeen 2J-NL cases, the length of the amplified regions varies from ten to several hundreds of kilobase pairs. This can be explained by the variable efficiency of gap repair synthesis, as observed in other assays in yeast and mammalian cell cultures. 32 The limited homology we observed at the final viral-chromosomal junctions (Supplementary Figure S2), is attributed to microhomologymediated primer extension events, as proposed in the MM-BIR (microhomology-mediated break-induced replication) model. 33 With similar molecular steps, the MJ cases may occur if the resolution of the viral integration events are not always efficient, allowing repetitive cycles of viral chromosomal insertions that lead to the formation of integrated concatemers, as represented in the MJ-CL category. An alternative, yet non-exclusive hypotheses, is that both ends of the invading virus behave independently. Subsequently, these ends may capture distant chromosomal sequences that ultimately create the complex set of junction events that are detected in the MJ-CL and MJ-SC categories. An important outcome of this scenario will be to copy and assemble together one or several physically distant chromosomal fragments. Thus, the identification of several host junctions in a single tumour would not necessarily mean that HPV integrated within several places and in several copies per genome, but rather integrated in a single locus that brought together a complex set of captured host sequences, upon multiple rounds of microhomology-mediated capture of chromosomal fragments. A recent study 20 found 3,666 HPV integration breakpoints in 103 cervical carcinomas or cell lines, but whether or not these breakpoints reflected single or multiple HPV insertions, remains to be established.
To estimate the impact of our HPV tumour classification system, we analysed the available patient clinical and biological data. The first unexpected finding was that, within the 2J-NL category, the majority (11/17) of the HPV18-associated cases were found in this group, whereas only 3/21 of the highly prevalent HPV16 genotype (here identified in 49% of all tumours) exhibited this signature ( Figure 4). Also, no HPV18 was observed in the EPI category ( Figure 4). In addition, the mean number of single-nucleotide variations within HPV16 is highest in the MJ-SC (17 nt) and EPI (14 nt) cases (Supplementary Figure S13A), with a different distribution of polymorphisms found among the groups along the entire sequence, as compared with the reference HPV16 sequence (gi|333031|lcl|HPV16REF.1) from the PaVE genome site (Supplementary Figure S13B). These polymorphic differences among the HPV16 cases in these groups may influence their oncogenic potential. The reason why the integration of HPV18 occurs according to a non-linear pattern associated with local amplification is not yet clear. Perhaps HPV16 and 18 genotypes have different mechanisms of integration (amplification versus deletion) due to their different regions and lengths of microhomologies between the viral and host genome. The second intriguing result concerned the striking age gap (up to 11 years) between the various HPV integrated and EPI categories (Supplementary Figure S14). The statistical comparison between these categories is not significant (Mann-Whitney U-test P value = 0.12), but since the presence of an episomal form of HPV is likely to be a prerequisite for integration, it is likely that the integrated versus episomal status is important, as a marker of oncogenicity levels that govern the length of time of neoplasias at their intraepithelial stages, before invasion develops. Oncogenicity could also influence the outcome of invasive carcinomas. We noticed in our case panel, that the rate of relapses was higher in tumours with a simple integration pattern, harbouring only two junctions than in cases with multiple viral insertions or with no insertions, developed in older women, an observation in line with the poor outcomes reported in tumours developed in younger women. 34 To test the robustness of the above mechanistic/clinical analyses, and comprehensively explore the significance of the multi-parametric genomic and functional diversity of the HPV virus, extended series of cervical HPVassociated cancers examined as reported here are warranted.
Altogether, our data sheds light on several key aspects of HPV-associated tumours from biopsy and ct-DNA samples, allowing new perspectives in clinical oncology. Most importantly, the optimal molecular characterisation of the associated HPV junction sequences by Capture-NGS will provide in 1-2 weeks the specific DNA biomarkers for the patient biological follow-up, during and after treatment, an approach that has recently demonstrated its interest in the follow-up of patients with breast cancers. 35 Collectively, the deep molecular analysis of viral-associated alterations nearby integration sites and genome-wide identification of candidate gene mutations should help the stratification of patients for innovative therapies in cervical cancer, targeting cellular genes 36,37 or viral sequences. 38,39 It should also be stressed that our approach allows the primary detection of any HPV-associated carcinoma using a simple blood test, an alternative method when detection of gynaecological cancers is difficult to implement. Another perspective will be to validate prospectively the potential value of this approach for the triage of cervical intraepithelial neoplasias at risk of progression. Finally, beyond cervical tumours, the present Capture-HPV approach is likely applicable to diagnose other HPV-associated cancers, such as head and neck carcinomas in which the proportion of HPV is increasing 40 and is also applicable for detecting other viral-associated tumours or the movement of mobile elements in genomes.

HPV probe design
Single-stranded biotinylated probes were designed with Roche NimbleGen (SeqCap EZ Developer targeting HPV, Roche NimbleGen Inc., Madison, WI, USA), in which over 22,000 individual and unique sequences, with an average probe length of 75 nucleotides, were synthesised to recognise 235 unique HPV types and variants, containing over 84 different low, medium, high-risk and even non-risk HPV's. The list of sequences used for the NimbleGen SeqCap EZ design of these biotinylated probes is shown in Supplementary Table S1. They were designed to recognise the entire 7.8-7.9 kb HPV episomal genomic sequences, allowing for all possible breakpoints to be detected in the virus. A masking strategy was used to reduce biases towards a non-uniform design containing many copies of highly similar sequences. Several additional unique chromosomal coding regions were added to serve as internal controls to estimate the viral load by measuring the ratio of HPV reads per kb against these similarly captured genes. The probes were additionally screened to contain no significant homologies to the human (hg19) genome.
Double capture of shotgun biopsy and circulating DNA libraries Two hundred to 500 ng of DNA were used for each patient sample for subsequent DNA library preparation. The DNA libraries are prepared by nebulising, which shears the DNA mechanically and randomly with nitrogen, polishing the double-stranded DNA ends and T/A cloning of adapters, containing unique barcodes for each sample. This is followed by the Double Capture system, using the SeqCap EZ Rapid Library Small Target Capture method, developed by Roche. This is adapted to enrich for and capture small DNA targets, in which the DNA library samples are hybridised, two successive times, with the target HPV oligonucleotide probes. The expected technical advantage of this double capture method versus a single capture 20 is to increase the specificity of the capture and enrich the fraction of relevant reads for subsequent bioinformatic analysis. DNA sequences captured by streptavidin were analysed using Roche 454 (GS Junior sequencing system) for patient cases 1-46, and Illumina MiSeq instruments, for the remaining 26 cases. Library samples were multiplexed, using one biotinylated HPV oligonucleotide probe set per 6 or 7 library samples. Single-end 454 shotgun reads were obtained, with an average length of about 500 nucleotides, allowing for easily identifiable viral-cellular junction reads, as well as purely viral reads, used for genotyping. The Illumina paired-end MiqSeq reads are 2 × 150 nt. Contig assemblies were also performed, using the Newbler software allowing for the full HPV integration status (the sequence of the viral insert with both 5′ and 3′ viral-cellular junctions). The sequencing parameters of each sample is provided in Supplementary Information.
The ct-DNA was prepared from 200 μl of serum collected from a patient blood sample at the time of diagnosis, using the QIAamp viral DNA (QIAGEN France SAS, Courtaboeuf, France) enrichment and blood kits. No sonication was performed before library preparation. The size profiles of the library DNA fragments were controlled on a Bioanalyzer before and after adaptor ligation. Fragments between 150 to 500 bps were purified and used for the subsequent Capture steps, performed as for the biopsy tumour samples described above.
Bioinformatic pipeline for analysis of viral-cellular DNA junctions and HPV genotypes Raw sequence data reads were obtained as either sff (binary) or fna (txt) formatted files. The sff files were converted into fastq files, using a python script: sff_extract 0.3.0 (http://bioinf.comav.upv.es/sff_extract/). All tools used thereafter were obtained from the Galaxy interface, at the Insitut Curie, whose development is based on a similar interface, available at (https://usegalaxy.org/) and is shown in Supplementary Figure S1.
We began by analysing the genotype of the sample. First, all of the sff reads were converted to fastq reads, that we aligned against a reference fasta file, containing 37 HPV low, medium and high-risk reference HPV genomes, obtained from the NIH PaVE website (http://pave.niaid.nih.gov/). Next, we used the Bowtie 2 tool from Galaxy set at default parameters, with an 'endto-end' mode of global alignment. The alignment (SAM file) is visualised using Tablet (http://bioinf.scri.ac.uk/tablet/). In this manner, we identify the HPV type of the sample, as Tablet displays the number of reads that are aligned against each of the 37 different HPV reference viral genomes.
Next, the viral-cellular DNA junctions were analysed. Here again we aligned the reads using the Bowtie 2 tool, against the HPV subtype reference genome, previously identified. However, we run a local alignment in order to obtain the purely HPV reads and also the hybrid reads containing at least 28 nucleotides aligned on HPV (the default value of Bowtie 2). We remove PCR-duplicates using the 'Mark Duplicate reads' Picard tool on the BAM file, to obtain unique, yet overlapping reads. We then convert the BAM into SAM files and filter all unaligned reads (using SAM tools).
The alignment with Tablet easily reveals simple (two clear junctions, 2J) from multiple viral insertion cases. For such 2J cases, the alignment of reads against the reference viral strain, allows the identification of non-pure HPV reads, or viral-host junctions, that do not entirely align against the viral reference genome. This alignment also shows the portion of the viral genome that is deleted upon integration. The CIGAR code is used to confirm the alignment ends of each NGS read at precise positions that contain an additional sequence, not aligned against HPV ('S' flag). In this manner, we obtain the non HPV portion of these reads that align against the hg19 reference using the blat tool (http://genome.ucsc.edu), to identify the precise insertion site of the virus in the human genome.
Alternatively, if a deleted or gapped region in the Tablet alignment is not detected, then we implement a filter to remove all purely viral reads and to highlight only viral-host junction reads. For this, we performed two local read alignments: the first against the HPV reference and the second against the hg19 reference sequence. After removing PCR duplicates and unaligned reads, the Galaxy tool, 'Compare Two Data Sets', is used and allows comparison between two data sets. Thus, by comparing both SAM files (HPV and Human) we find the common read names that are partially aligned on both HPV and human sequences, that we call hybrid reads. As before, we visualise the SAM files in Tablet, with clearly visible viral-host breakpoints.

Assignment of viral-cellular sequences and breakpoint coordinates
The coordinates of the HPV-genomic breakpoints were assigned using the hg19 reference genome and HPV Pave website, for reference HPV genomes. Fifty-to 55-nt sequences overlapping the viral-cellular breakpoints were aligned, resulting in 5 up to 24 nts of gapped sequence identity either at or surrounding the breakpoint, within the entirely aligned sequence (Supplementary Figure S2). Exact overlapping sequences at the breakpoint, were assigned to the human genome.

De novo read assembly (Newbler)
We performed a de novo assembly, when possible, using the Newbler software from Roche. The raw data (sff of fastq file) was used as input and the following default parameters were used as following: minimum read length 20; minimum overlap length 40; minimum overlap identity; 90%.
The output fna file, containing all of the contig sequences were aligned against the HPV genome. Then, the non-HPV portion of the contig was identified by aligning on hg19 (with the BLAT tool). This allowed us to confirm junctions that we identified previously, from single reads and confirm our integration breakpoints obtained from our single reads. Indeed, for most of the 2J cases, Newbler was able to assemble one or two long contigs, encompassing both the 5′ and 3′ flanking viral-human breakpoints, as well as the entire integrated viral sequence.

PCR validation of viral-host junctions
The DIPS-PCR (Detection of integrated papillomavirus sequences) method, followed by Sanger sequencing was performed on all samples, to validate our NGS identified junctions using the Capture-HPV-NGS method, as described previously. 25,26 Measurement of viral loads HPV viral loads were determined in two ways. First, as previously described, 8 by quantitative PCR, in which the HPV E7 gene is quantified relative to the reference PSA gene (2 copies of DNA per cell) using quantitative PCR (Sybr Green PCR Core Reagents). The viral charge is obtained using the formula CV = 2 −ΔCt × 2. Differently, the HPV viral load was measured from the capture NGS results, by counting and comparing the nucleotide coverage of the HPV sequences and those of the three control house-keeping genes, included in the same tube as the HPV capture probes. These three unique genes (GAPDH, KLK3 and RAB7A) with a size (6-11 kb each) comparable to the HPV were captured with a total of~3,000 copy-number control probes. The estimated viral load was determined by the ratio of HPV versus control gene read coverage, normalised per nucleotide. Correlation between these methods is shown in Supplementary Information.

Genome-wide analysis of HPV integration using CytoScan Affymetrix
All samples (except the #45 and 46 cases, with no detectable HPV DNA) were analysed at the Genetic platform at the Institut Curie, using the Genome-Wide Human CytoScan HD Array from Affymetrix, for comprehensive whole-genome coverage of DNA from cervical biopsy samples (includes 2.67 million markers for copy-number analysis, with 750,000 SNP probes and 1.9 million non-polymorphic probes). Analysis was performed using the Chromosome Analysis Suite software (ChAS) version 3.0 (Affymetrix), as well as the Genome Alteration Print (GAP) software, developed at the Institut Curie. 41 Frequency of Amplicon, Gain and Loss (FrAGL) was performed on the bioinformatics platform (U900, INSERM) of the Institut Curie.