Blood samples from 10 patients treated under a systematic surgical algorithm at The University of Texas MD Anderson Cancer Center (Houston, TX) were retrieved from the Gynecologic Tumor Bank after written informed consent was obtained under a protocol approved by the Institutional Review Board (LAB10-0850), as we described previously (20, 21). These patients had upfront inoperable disease and received carboplatin-based NACT (neoadjuvant chemotherapy). After 3-4 cycles of NACT, patients were considered to have excellent response (NACT-ER) if they had a complete response or only microscopic disease at the time of interval surgery; they were considered to have poor response (NACT-PR) if they had stable or progressive disease on radiologic evaluation and/or suboptimal interval cytoreduction after NACT, according to Response Evaluation Criteria in Solid Tumors version 1.1.

Genomic DNA was extracted with the QIAamp DNA Mini Kit (Qiagen) and then quantified using the Quant-iT PicoGreen dsDNA Assay Kit (ThermoFisher Scientific); quality was assessed using Agilent High Sensitivity D5000 ScreenTape and Reagents on the TapeStation 4200 system (Agilent Technologies). cfDNA was extracted using the MiniMax High Efficiency Cell-Free DNA Isolation Kit (#A17622-50; Apostle) and then quantified using the Quant-iT PicoGreen dsDNA Assay Kit; quality was assessed using Genomic DNA ScreenTape and Reagents on the TapeStation 4200 system. Each genomic DNA sample (up to 200 ng, based on the PicoGreen quantification) was sheared (mechanically fragmented) using the E220 focused-ultrasonicator (Covaris) with the following settings: peak incident power, 200 Watts; duty cycle, 25%; cycles per burst, 50; duration, 10 seconds; iterations, 70. To ensure the proper fragment size, we examined the samples on the TapeStation 4200 system using Agilent High Sensitivity D1000 ScreenTape and Reagents. Libraries of the sheared DNA and unsheared cfDNA were prepared using the SureSelect XT_HS2 DNA Reagent Kit (Agilent) with 384 unique dual sample indexing and dual molecular barcodes to better suppress false positives and more accurately detect low VAFs. This protocol consisted of 3 enzymatic reactions for end repair, A-tailing, and adaptor ligation followed by barcode insertion by polymerase chain reaction (PCR) using Herculase II Fusion DNA Polymerase (Agilent; 8-14 cycles, based on input DNA quality and quantity). PCR primers were removed by using 1x volume of the AMPure PCR Purification kit (Agencourt Bioscience Corporation). The quality and quantity of the prepared libraries were evaluated using the High Sensitivity D1000 ScreenTape on Tapestation 4200 system to verify correct fragment size and ensure the complete removal of primer dimers. Subsequently, the prepared libraries were individually hybridized to SureSelect Human All Exon V4 probes (Agilent). The hybridization steps were automated on the Sciclone G3 NGSx Workstation (PerkinElmer, Inc.). Agilent Captured regions of interest were hybridized as single-sample reactions using 500–1000 ng of the prepared library as input. All hybridization and post-hybridization captures and washes were performed according to Agilent’s protocol. Briefly, the capture reagents and probes were added to the prepared libraries, and the mixture was incubated at 65°C in a thermocycler with a heated lid for up to 24 hours. The targeted regions were captured using streptavidin beads, the streptavidin-biotin-probe-target complex was washed, and the captured libraries were enriched by PCR amplification according to the Agilent’s protocol. The quality and quantity of each captured sample were analyzed on the TapeStation 4200 system using the DNA High Sensitivity Kit. The captured libraries were sequenced on the Illumina NovaSeq 6000 platform for 2 × 150 paired-end reads with an 8-nt read for indexes using Cycle Sequencing v1.5 reagents (Illumina). Demultiplexing was performed using Illumina’s bcl2fastq or BCL Convert software to generate paired-end reads based on the dual indexes and remove sequences with incorrectly paired P5 and P7 indexes. The Agilent Genomics NextGen Toolkit (AGeNT) was used for molecular barcode extraction and trimming.

Sequencing output fastq files were first processed with fastp (version 0.23.0) using parameters –disable_adapter_trimming, –qualified_quality_phred 20, –unqualified_percent_limit 30, and –average_qual 20. Processed reads were then mapped to human genome hg38 using the Burrows-Wheeler Aligner (version 0.7.15) (22) using command “bwa mem -Ma” and default settings. The result file is converted to BAM file using samtools (v1.9), and genome coverages were estimated using “samtools idxstats” to obtain an average coverage on each chromosome followed by calculating a median weighted coverage adjusted by chromosome sizes as reported in Agilent_Human_Exon_V4_Regions.S043800110.hg38 bed file. To minimize PCR introduced sequencing error, we edited the BAM files and collapsed reads with the same molecular barcodes. When more than two sequencing reads are mapped to a same chromosomal start and end position, they are collapsed into two reads with the same sequences which is the most occurred sequence at this position. Only reads with properly positioned molecular barcodes were handled (ca. 85%). Reads with improper barcodes or mapping were discarded. In addition, paired reads with overlapping sequencing area are clipped on the forward strand to avoid counting the same sequence twice. The duplication percentages of original and UMI collapsed BAM files were reported by Picard MarkDuplicate function (Picard v2.27.4) (21). Varscan (version 2.4.2) was used to call variant SNPs on a pile up file generated by samtools with parameters –min-coverage 20, –min-reads2 2, –min-avg-qual 20, –p-value 0.05, and –min-var-freq 0.005. Result variants were annotated using GATK Funcotator (v4.2.4.0, data source v1.7.20200521g) and default settings. CHIP variants were defined as heterozygous variants found in both the WBC and cfDNA samples from the same patients with VAF of less than 5% in the WBC sample and more than 0.5% in cfDNA samples. Variants not meeting these criteria were considered non-CHIP variants. Known CHIP genes (n = 50) were compiled from previous studies (2, 23). Cancer Gene Census genes (n = 733) were collected from COSMIC (https://cancer.sanger.ac.uk/census) (24). Hypergeometric tests were used to examine whether Cancer Gene Census were enriched within the CHIP and non-CHIP gene sets.

Comparative analysis of the differential CHIP variants among clinical groups (all pre-NACT samples versus all post-NACT samples; PR pre-NACT samples versus PR post-NACT samples; ER pre-NACT samples versus ER post-NACT samples) were performed using R limma package (R v4.1.0, limma v3.52.4). We considered P values less than 0.05 to indicate statistical significance. Owing to the limited number of samples, most of the adjusted p-values obtained from limma analysis did not indicate significance. We therefore used an additional Wilcoxon test to prioritize limma results for putative top differential CHIP variants. For pathway analysis IPA (ingenuity pathway analysis) from Qiagen was used.

The study included 9 patients reported previously (20, 21). Of these patients, 4 had an excellent response (ER) to NACT (NACT-ER), and 5 had a poor response (PR) to NACT (NACT-PR) ( Table 1 ). Isolated cfDNA, white blood cell (WBC) DNA, and tumor DNA from each patient were subjected to ultra-high-depth whole-exome sequencing using unique molecular barcodes ( Figure 1 ).

We defined low-frequency candidate CHIP variants as heterozygous variants that were present in both WBC DNA and cfDNA and that had a VAF below 5% in the WBC DNA and above 0.5% in the cfDNA. We considered genes in which CHIP mutations were previously identified to be known CHIP genes (2, 23). The total number of low-frequency CHIP variants was 93,088, involving 13,780 genes ( Supplementary Table 1 ). Of these variants, 463 were in 44 known CHIP genes, and 5,548 were in 598 genes listed in the Cancer Gene Census in the Catalogue of Somatic Mutations in Cancer (COSMIC) ( Supplementary Table 2A ). The total number of recurrent candidate CHIP variants present in at least 2 patients was 47,961, involving 13,509 genes. Of these variants, 266 were in 44 known CHIP genes, and 2,874 were in 585 genes listed in COSMIC ( Supplementary Table 2B ).

After quantifying the candidate CHIP variants in the liquid biopsy samples, we compared all patients’ pre- and post-chemotherapy samples. Post-chemotherapy samples consistently had fewer CHIP variants than pre-chemotherapy samples did ( Figures 2A, B ). To determine whether this change depended on the variant’s original allele frequency (before chemotherapy), we segregated the CHIP variants into 4 groups based on their VAFs in cfDNA in pre-chemotherapy samples (<1%, 1-2%, 2-5%, and >5%) and compared the pre- and post-chemotherapy frequencies of the CHIP variants in each group. This difference in number across groups, in terms of decrease in post-chemotherapy samples, was significant (p= 0.002, Figure 2B ). However, after segregation, the total number of CHIP mutations whose VAFs were between 2 and 5% did not differ significantly between pre-NACT cfDNA and post-NACT cfDNA ( Figures 2C–F ).

To identify the candidate CHIP variants that were potentially selected by chemotherapy, we performed a comparative analysis of all patients’ pre- and post-chemotherapy samples to identify variants whose VAFs increased after chemotherapy. For this analysis, we first considered each CHIP variant’s VAF to be its average VAF in both WBC DNA and cfDNA. We identified 77 mutations with increased allele frequencies after chemotherapy ( Figure 2G , Supplementary Table 3A ), but none of these increases were statistically significant after adjustment for multiple hypothesis testing. We then sought to identify the pathways that were enriched among the genes whose CHIP variants had higher VAFs after chemotherapy. We identified enrichment in the activation of NMDA receptors and postsynaptic events (p = 8.93E-05; 4.7% genes overlap), transcriptional regulation by RUNX3 (p = 1.43E-04; 4.2% genes overlap), ERK/MAPK signaling (p = 3.52E-04; 2.3% genes overlap), Signaling by TGFBR3 (p = 3.54E-04; 6% genes overlap), and transcriptional regulation by MECP2 (p =7.00E-04; 4.8% genes overlap). In addition, the genes whose CHIP mutations had higher VAFs after NACT were involved in cellular and molecular functions, including cellular development (p = 1.36E-02-9.64E-06), cellular growth and proliferation (p = 1.36E-02-.5 = 9.64E-06), Cell Death and Survival (p = 1.31E-02 4.6E-05- 0.00329), gene expression (p = 1.36E-0.2-5.99E-05), and Cell Morphology (p =1.24E-02- 9.52E-05). The candidate upstream regulators for some of these genes were TNKS1BP1 (p= 1.78E-05), DNMT3A (p = 5.95E-05) and TET3(p = 2.51E-04), some genes previously associated with leukemia (25, 26). The top candidate CHIP mutations included a missense mutation (p.R708L) in CDK12, which encodes cyclin-dependent kinase 12, a protein involved in the tumorigenesis of several cancers. To test the reliability of circulating DNA as a source material for CHIP variant quantification, we performed the same comparative analysis using only the candidate CHIP variants’ VAFs in cfDNA; of 78 variants with higher VAFs after chemotherapy, 52 were also identified in the previous analysis ( Supplementary Table 3B , Supplementary Figure 1A ).

We then sought to determine whether CHIP variants’ VAFs differed between pre- and post-NACT liquid biopsy samples for patients with NACT-ER and for patients with NACT-PR. For PR patients, 85 CHIP mutations had higher VAFs in post-NACT samples than in pre-NACT samples, but none of the differences was statistically significant ( Figure 3A , Supplementary Table 4A ). For ER patients, 2 candidate CHIP mutations had significantly higher VAFs after NACT: p.V100G, in PRKAR1B (adj. p = 0.001845171); p.E92A, in DUSP9 (adj. p = 0.008174934) ( Figure 3B , Supplementary Table 5A ). We then compared CHIP variants’ VAFs in cfDNA before and after chemotherapy. In PR samples, 79 variants with higher VAFs after chemotherapy ( Supplementary Table 4B , Supplementary Figure 1B ). In ER samples 75 variants had higher allele frequency after chemotherapy ( Supplementary Table 5B , Supplementary Figure 1C ).

A comparison of CHIP variants’ VAFs in post-NACT samples from NACT-PR and NACT-ER patients revealed that 59 variants had higher VAFs in PR samples than in ER samples ( Supplementary Table 6A , Supplementary Figure 2A ). A comparison of CHIP variants’ allele frequencies in post-NACT cfDNA samples from ER and PR patients revealed that 45 variants had higher allele frequencies in PR samples; of these, 28 were also identified in the previous analysis, but these differences were not statistically significant ( Supplementary Table 6B , Supplementary Figure 2B ).

To identify novel prognostic markers that can be detected with non-invasive techniques (such as liquid biopsies) in ovarian cancer patients, we sought to identify candidate CHIP variants whose frequencies were higher in pre-NACT samples from NACT-PR patients than those from NACT-ER patients ( Figure 3C ). Our analysis identified 1 CHIP variant, p.Y73*, in MRPL50 (adj. p = 0.018876) ( Supplementary Table 7A ). MRPL50 encodes mitochondrial ribosomal protein L50, and the detected variant is not a previously known CHIP variant or in the COSMIC database. A comparison of CHIP variants’ allele frequencies in pre-NACT cfDNA samples from PR and ER patients revealed 70 variants whose allele frequencies were higher in PR samples ( Supplementary Table 7B , Supplementary Figure 2C ).

To investigate the intra-tumoral presence of clonal hematopoietic cells characterized by CHIP mutations, we used deep sequencing to analyze DNA from pre- and post-NACT tumor samples matched to the cfDNA and WBC DNA. We found that 26-46% of the CHIP mutations in cfDNA matched those in tumor DNA ( Figure 4A , Supplementary Table 8 ). The LIMMA test identified several candidate CHIP mutations whose VAFs were higher in post-NACT tumor DNA than in pre-NACT tumor DNA, but these differences were not significant after p-value adjustment ( Figure 4B , Supplementary Table 9 ). Among patients with NACT-PR, 221 CHIP variants had higher allele frequencies, and 93 had lower allele frequencies, after chemotherapy ( Supplementary Tables 10A, B ). Among patients with NACT-ER, 165 CHIP variants had higher allele frequencies, and 121 had lower allele frequencies, after chemotherapy ( Supplementary Tables 11A, B ). A pathway analysis of the genes affected by variants with higher allele frequencies after chemotherapy identified the Myelination Signaling Pathway (p = 1.24E-04; overlap 2.7%) as the top canonical pathway.