DSB-TRIP works by random integration of a specially designed DSB repair pathway reporter into hundreds or thousands of genomic locations in a pool of cells, by means of a transposon vector. The reporter is short [∼650 base pairs (bp)] and devoid of transcriptionally active sequences that could change the local chromatin environment (Figure 1A). Each copy of the reporter is marked by a random barcode, which allows for decoding of individual reporters and linking them to their genomic location. First, the genomic locations of the reporter integrations in the pool of cells are mapped. It is also possible to generate clonal cell lines that carry up to dozens of reporters. Next, a DSB is introduced inside each reporter by means of Cas9. Repair of the resulting DSBs results in specific “scars” [insertions and deletions (indels)] that can be used to identify the repair pathways that were active at a DSB (Allen et al., 2018; Chakrabarti et al., 2019; Chen et al., 2019; Shen et al., 2018; van Overbeek et al., 2016; Brinkman et al., 2018). Our current reporter can detect NHEJ, MMEJ and SSTR. After DSB induction and DNA repair, the genomic DNA is extracted and the reporter “scar” and the flanking barcode are jointly amplified by PCR and subjected to high throughput sequencing. A computational pipeline counts the scars for each individual barcoded reporter, infers from these counts the relative activity of each pathway, and links the results to the genomic location. Overlay with epigenomic mapping data then uncovers any correlations between pathway usage and local chromatin features.

A key feature of DSB-TRIP is that an identical reporter sequence (except for the short barcodes) is integrated into many different chromatin environments. Hence, differences in pathway usage between integration sites can be attributed to differences in the chromatin context. Thus, insights are obtained in the impact of the chromatin environment on DSB repair.

The protocol consists of two main components: the wet-lab experiments, and the computational analysis of the resulting sequencing data.

The wet-lab part starts with the design of the DSB repair reporter, which is a DNA sequence that, when repaired after a DSB, produces a specific indel pattern that can be associated with either NHEJ or MMEJ. We then describe how to make reporter plasmid libraries, how to transfect them into K562 cells and how to generate cell pools and clones carrying multiple integrations (Figure 1B). Next, two sets of data are generated: mapped reporter integrations and indel scores per reporter. First, the genomic location of the integrated pathway reporters (IPRs) in these clones are mapped by inverse PCR (iPCR), as explained in detail in (Akhtar et al., 2014). Second, the DSBs are induced by Cas9 induction, cells are cultured up to 72 h to allow the DSBs to be repaired, after which genomic DNA is collected to create sequencing libraries.

The sequencing data of the mapping and the indels are then processed with a computational pipeline that produces tables of the mapping per barcode and indels per barcode for further analysis.

The most basic implementation of the data analysis pipeline uses two types of sequencing data: iPCR reads used to link barcodes to the integration locations in the genome (mapping) and indel PCR sequences to assess repair outcomes. The standard output for the iPCR is a table containing the barcode, its genomic location (chromosome, start, end, and orientation), reads mapped, mapping quality and mapped sequence. The output for the indel PCR is a table containing four columns: the barcode sequence, the mutation class (e.g., wild-type, insertion, deletion or “unclear”), the size of the detected indel in bp and the number of occurrences. Two main extensions can be implemented for more detailed output. First, the pipeline can be run to count specific repair scars. This can be used to differentiate MMEJ from non-MMEJ scars with the same indel size. This functionality is available in the pipeline but requires a significantly longer runtime. The second option available, is the ability to add a recognition sequence to detect homology directed repair.

A prerequisite for DSB-TRIP is a cell line in which a DSB can be induced in a specific and reproducible way. We recommend establishing a clonal founder cell line with a stably integrated, inducible Cas9. We use K562 cells expressing Cas9 from the human PGK promoter and fused to a destabilization domain (DD) at its N terminus. This DD causes active degradation of the protein unless it is stabilized by addition of the small molecule Shield-I (Banaszynski et al., 2006; Brinkman et al., 2018). This allows for tight control of Cas9 activity. We note that in our cell line (clone K562#17) full induction of DD-Cas9 takes ∼16 h, which should be taken into account in the experimental design. Possibly, Cas9 can be introduced into cells by lentiviral transduction, or by transient transfection with an expression vector or Cas9 ribonucleoprotein. However, this may not yield homogeneous expression in 100% of the cells, which may compromise the data quality.

Instead of K562 cells, other cells may be used. However, there are some practical restrictions. First, the cell line must tolerate DNA transfections and transposon integrations; this is a requirement for the reporter insertions as well as for sgRNA transfections. Second, the cells must tolerate some levels of DNA damage, as some individual cells in the pool might carry more than 20 IPRs in their genome, and thus may potentially need to cope with as many simultaneous DSBs. Third, epigenome mapping data must be available for the cell line. This is essential for linking of the DSB-TRIP results to the chromatin state of the IPRs.

One should consider the following elements when designing the reporter. First, a reporter of small size (≤1 kb) and devoid of any active transcription will have less impact on the chromatin state at the integration site. Second, the reporter sequence should produce signature “scars” after repair of the DSB that are characteristic of the repair pathway(s) of interest. We used a short sequence derived from the LBR gene (Schep et al., 2021), but it may be replaced by other sequences depending on the pathway of interest. Third, when designing the reporter sequence it is important to have the correct restriction sites at each end of the sequence (Section 1) as well as the IndelPCR1_rv primer binding site at the NheI end to be able to follow the library preparations steps (Section 13).

Despite the power of DSB-TRIP, it has some limitations. First, DSB-TRIP is based on the detection of signature indels that are produced by the respective pathways. However, HR does not usually generate indels, and hence this pathway cannot be measured with DSB-TRIP. The reporter that we present here is also sensitive to large deletions or extensive resection that might remove primer binding sites or the reporter barcode. In part this can be overcome by using Tn5 transposon-based library preparation as we presented in (Schep et al., 2021) or with UdiTaS (Giannoukos et al., 2018) which requires only one specific primer site. Second, while Cas9 has been the most prominent tool for gene editing, the DSBs that it creates may not be representative for DSBs that occur naturally. For example, the cutting and repair rates that have been observed with Cas9 appear to be extremely slow (Brinkman et al., 2018). This should be kept in mind while interpreting the data. Due to the slow rates of cutting and repair (a typical DSB-TRIP experiment takes 3 days), it is difficult to link observed pathway activities to cell cycle stages. Finally, in the cell pools the number of integrations per cell can vary substantially, which may lead to different DNA damage loads that in turn may affect the DNA repair kinetics or pathway balances.

While the generation of TRIP cell pools and clones can be time consuming there are many ways to use these TRIP cell pools and clones. The cell pools offer many different genomic loci, but note these experiments also require many cells to have sufficient coverage per IPR. For more multiplexed applications, such as small molecule or CRISPR screens, we recommend the use of clones carrying multiple IPRs. DBS-TRIP is compatible with experiments in 96- or even 384-well plates. This format is perfect for screens and other automated applications. For instance, a multiplexed, automated time series analysis in 96-well format was presented in Schep et al. (2021).

The reporter assay can be modified in many ways. For example, the DSB target site may modified with different microhomologies to study MMEJ in more detail. It is also possible to add variations in the target sequence to study how chromatin might affect slight changes in affinity of Cas9. This may provide a better understanding of off-target cutting by Cas9 and the ensuing repair (Figure 2A, middle). Using the exact same reporter as presented here it is feasible to study SSTR in the context of chromatin, by co-transfecting a single-stranded oligonucleotide that carries a small insertion (Schep et al., 2021). It may also be possible to provide different types of templates, for instance a double-stranded DNA template (Wienert et al., 2020). Furthermore, the effects of adding mismatches and varying the homology arm lengths may be investigated (Figure 2A, bottom) (Richardson et al., 2018; Hussmann et al., 2021). Potentially, time-series measurements on DSB-TRIP cell pools or clones combined with mathematical modelling (Brinkman et al., 2018) may reveal how rate constants of individual steps of the cut-and-repair process are affected by the local chromatin state.

It will be interesting to apply DSB-TRIP with other CRISPR endonucleases such as Cas12a (Zetsche et al., 2015) or rare-cutter restriction enzymes such as I-SceI (Niu et al., 2008) or AsiSI (Iacovoni et al., 2010) (Figure 2A, bottom). I-SceI was successfully applied to study Single Strand Annealing in a TRIP assay (Pokusaeva et al., 2021). Another option are transcription activator-like effector nucleases (TALENs) (Boch et al., 2009; Gaj et al., 2013; Moscou and Bogdanove, 2009), especially since it was found that these nucleases are more efficient in heterochromatin, compared to Cas9 (Jain et al., 2021). Furthermore, mismatch repair may be studied using base editors (Figure 2A, top) [reviewed in (Anzalone et al., 2020)]. For each of these endonucleases it is essential that they are expressed in a tightly controlled inducible manner, otherwise the DSBs are already generated and repaired while the cell pools with IPRs are being established.

Finally, it will be interesting to study chromatin context effects on DNA repair in cells with particular genetic alterations, such as mutations in specific repair proteins or defects in chromatin organisation, such as Hutchinson-Gilford Progeria that is the result of specific mutations in the Lamin A gene (De Sandre-Giovannoli et al., 2002; Eriksson et al., 2003). With its flexibility, multiplexing ability and detailed readout, DSB adds a versatile tool to study DNA repair, gene editing and the interplay between repair pathways and chromatin.

The DSB-TRIP pipeline is written as a Snakemake pipeline (Molder et al., 2021). Two files are required for successful execution: 1) a meta-data file with information on the samples, and 2) a pipeline configuration file with pipeline parameters.

1. Meta-data file

A tab-delimited file with information on the samples: sample name, file location, type of PCR (e.g., iPCR, indelPCR), and plasmid guide information.

2. Configuration file

The configuration file determines the DSB-TRIP pipeline execution. Parameters include but are not limited to: description of the read structure, thresholds used in the pipeline, information on the reporter for DNA scar analysis, and alignment parameters.

The DSB-TRIP pipeline consists of several modules which will perform the following four tasks: pre-processing of the sequencing reads, barcode clustering, genomic alignment of the iPCR reads and scar analysis of the indel PCR reads.

1. Barcode retrieval

First, the barcode is extracted from iPCR and the indel reads. This implementation is based on Cutadapt (version 1.11) (Martin, 2011) and custom scripts. For both read types, a 16 bps barcode is positioned within a constant 5′ adapter sequence: “GTCACAAGGGCCGGCCACAA{barcode}TGATC”. The barcode is extracted using a three-step strategy: 1) matching of the entire sequence, 2) matching of the 5′ sequence, 3) matching of the 3’ sequence. This approach ensures accurate retrieval of barcodes shorter or longer than 16 bps (15 or 17 bp barcodes can be present). Only reads with successful barcode retrieval are selected for further processing. This barcode retrieval is performed on read 1 of the iPCR paired-end sequencing reads. An additional constant sequence is required for read 2 of the pair. Finally, the reverse complement of these constant patterns are removed from the 3′ end of the reads if present due to read-through into the backbone.

2. Barcode clustering

The retrieved barcode sequences are classified as “genuine” or “mutated” using the starcode clustering algorithm (Zorita et al., 2015). Barcodes are clustered based on Levenshtein distance (<2) and the most abundant barcode in the cluster is considered the “genuine” barcode. Currently this step is used as a filtering step, but alternatively, barcodes can be rescued by assigning all information from the mutated barcodes to the genuine barcode.

3. Reporter alignment

Remaining iPCR sequences are aligned to the human genome with Bowtie 2 (Langmead and Salzberg, 2012). This is carefully done to deal with potential chimeric reads created by the circular ligation [as shown in (Akhtar et al., 2014) Supplementary Figure 1]. For such reads, two pieces of genomic DNA are ligated at an unknown junction. Briefly, to remove these events, two sequential alignment steps are performed. First, the read pairs are mapped independently and the mates reading directly into the genome from the edge of the transposon are filtered based on potential soft-clip at the 5’ end. These reads are either discarded, or if >16 bps soft-clipping occurred, the clipped sequence is extracted and realigned to the genome in a second alignment step. This ensures that in case of chimeric reads, the genomic sequence found next to the integrated reporter is aligned instead of some distant genomic sequence that got inserted at the circular ligation step. Reads that pass the quality thresholds are aggregated based on their barcode sequence. For every barcode, read counts and mapping quality for the two locations with most supporting reads are returned. This allows custom filtering to select trustworthy integration sites.

4. DNA scar analysis

The DNA scar analysis tests for two points: the presence of the wild-type sequence and the size of a possible indel compared to the original sequence. First, reads are matched against known recognition sequence(s), made up of the wild-type sequence and optionally, any number of expected introduced mutations. Second, the location of specific downstream sequences is compared with their location in the wild-type fragment to calculate an indel size. In our case, these are seven 6 bps sequences that are originally are positioned between 83 and 124 bps from the start of the sequence (Figure 3D, iii). This approach allows us to pick up the complete range of deletion sizes that can be captured by the PCR performed and deal with potential read errors or mutations downstream of the break site.

When a recognition sequence is found, the call will be named after the recognition sequence (“wt” for wild-type sequence, or in case of a directed mutation, the name given in the configuration file). When no recognition sequence is found, the call depends on the indel size, this call will be “del”, “wt_point_mut” or “ins” for indel sizes <0, 0, >0, respectively. When no recognition sequence or indel size can be calculated, the call will be “not_clear”. The indel size will be calculated and reported even if in the first step a recognition sequence was found. This way, some reads might be given a “wt” call, but an indel size other than 0 (e.g., −1). This can occur due to either a technical artefact or biological mutation outside of the Cas9 targeted mutagenic event of interest.

Two alternatives to this scar analysis are implemented: 1) pairwise alignments that capture the exact DNA scar, 2) direct comparisons of the DNA scars with expected scars as predicted by inDelphi or FORECasT (https://partslab.sanger.ac.uk/FORECasT) (Allen et al., 2018; Shen et al., 2018). Both of these options will allow more fine separation of repair scars at the cost of significantly longer processing time. This allows for more detailed analysis such as separating MMEJ scars from NHEJ scars that resulted in the same size deletion.

The DSB-TRIP pipeline creates various output files and log files with details on pipeline execution. Output is organized following the four steps described previously.

1. Barcode retrieval

The statistics file for every sequencing sample is produced and includes for each read all the intermediate steps in the barcode retrieval analysis. This information is useful to identify and solve read parsing issues.

2. Barcode clustering

The barcode clustering returns a table with three columns containing barcodes: one column with clustered barcodes, the second with the most prominent within the cluster (we call this the “genuine” barcode), and one column with the total number of occurrences of barcodes within the cluster. The two columns containing barcodes are used to filter intermediary files for iPCR and indel reads that only the rows with “genuine” barcodes are left. The total cluster count is ignored by the pipeline, but can be used to assess how well this barcode clustering performed.

3. Reporter alignment

For each iPCR sample two tables are generated (we call these the “mapping tables”) containing barcode locations and quality measures, one for each mate of the paired end reads. The P7 read (read 2, file ending in “2.table”) is most important since that read contains the part of the sequence starting at the end of the transposon (Figure 3D, ii). The P5 read (read 1, file ending in “1.table”) is partly made up of the barcode and surrounding transposon sequence, there is still some genomic DNA left in the sequence which should be aligned at a nearby GATC side where DpnII cut the DNA before circularization. The table of this mate in the read pair can be used as quality control.

4. DNA scar analysis

The “.count” files generated by the scar analysis contain a four column table with the total number of reads per combination of barcode, mutation call and indel size produced by step 4 of the pipeline (see above).

Finally, the various pipeline output files can be used to evaluate the DSB-TRIP experiment. We advise several quality controls and downstream filtering steps, as described below:

The first thing to do is to select the barcodes which give an accurate and unique location. We advise to select all barcodes which are supported by at least 5 iPCR reads with an average mapping quality larger than 10 at the primary location, having at least 95% of the reads located at this locus, with not more than 2.5% of the reads at a secondary location. These ratio’s and averages can be calculated from the read counts and sum of mapping quality found in the “mapping table” (file ending in “.2.table”, see above). This will remove barcodes with multiple alignment locations either because of technical reasons such as repetitive sequences or because a single barcode is integrated in multiple locations.

After quality control and filtering, indel ratios can be calculated. In (Schep et al., 2021) this was done by normalizing each replicate over library size and biological replicates were averaged. The frequency of each indel type as proportion of total reads was calculated on that average. Pathway frequency per IPR was calculated as a proportion of the specific mutation over all indels (excluding wild-type sequences).

It is important to consider that wild-type (wt) sequences can be of various origins during the analysis. They can be either uncut or perfectly repaired by HR, NHEJ or another pathway. Additionally, wt sequences could undergo cycles of cutting and repair as long as they are repaired perfectly. However, we previously estimated that perfect repair of this reporter sequence is rare (Brinkman et al., 2018). As many factors can affect the abundance of the wt sequence in the cells, it is crucial to not overinterpret these results and take cell viability, transfection efficiency and various repair pathway activities into account.

The next step is to classify indels in pathways. Schep et al., 2021 used multiple gRNA target sites within the LBR gene, the main one was termed “LBR2” (“LBR1” in Brinkman et al., 2014) (Brinkman et al., 2014; Brinkman et al., 2018; Schep et al., 2021). This target site was used since it had a predictable repair outcome with very prominent outcomes for MMEJ events (7 bp deletions, flanked by 3 bp microhomologies) and NHEJ events (1 bp insertions). 7 bp deletions were classified as MMEJ, and 1 bp insertions were classified as NHEJ. Schep et al., 2021 focused on the balance between these outcomes (Brinkman et al., 2014; Brinkman et al., 2018; Schep et al., 2021). This classification can easily be adapted for different guides, when clear microhomologies and NHEJ outcomes are known. To discover which indels are caused by the NHEJ pathway, in Schep et al., 2021 the M3814 inhibitor was used (Schep et al., 2021). Indels with an indel ratio of at least 0.01 in the DMSO setting and a significant decrease (adjusted p-value < 0.05) in the M3814 inhibitor setting (one-sided Wilcoxon test) were classified as NHEJ. Alternatively, for guides that result in a broader spectrum of multiple mutations, a pipeline setting can be used to detect specific MMEJ events. This alternative pipeline setting performs pairwise alignments with the wild-type sequence was used, resulting in more detailed scar information at the cost of extra computational time. Adapting this strategy for LBR2, other microhomologies can be identified. These include a 14 bp deletion with a 3 bp microhomology and a 22 bp deletion with a 6 bp microhomology.

After mutation counts are separated by MMEJ and NHEJ (and unclassified) events, different calculations can be made for pathway balance and mutational rates. To calculate basic mutational rates a division over total indel can be a simple measure. For comparison in pathway balance between NHEJ and MMEJ, a relative pathway proportion can be calculated using: M M E J ( M M E J + N H E J ) .