Establishment and genomic characterizations of patient-derived esophageal squamous cell carcinoma xenograft models using biopsies for treatment optimization

Patients and tissue samples

The research proposal had been approved by the Medical Ethics Committee of Peking University Cancer Hospital according to principles of the Declaration of Helsinki. Written informed consent was obtained from all study participants for their information to be stored in the hospital database and used for future research at the time of follow-up ascertainment. All the patients with pathologically confirmed esophageal carcinoma and available endoscopic biopsied samples were included in the study.

Establishment of PDX models

NOD/SCID mice (6 weeks) were from Beijing Vital River Laboratory Animal Technology Co., Ltd. Two tissue fragments (~ 2 × 2 × 2 mm³/fragment) from fresh biopsies were obtained from each patient (P0 = passage zero), and were subcutaneously implanted into the flank of each mouse under sterile conditions. Tumor growth was assessed by palpation or Vernier calipers twice weekly. The established PDX model was passage 1 (P1). Mice were euthanized and harvested fresh tumor fragments were re-implanted into other mice when P1 tumors reached ~ 750 mm³. If animals showed disease, the tumor was collected when palpable. Subsequent passages were P2, P3, and P4. All procedures were performed under sterile conditions at Peking University Cancer Hospital specified-pathogens free facility and carried out in accordance with the Guide for the Care and Use of Laboratory Animals of the NIH. Each model derived from individual patients was passaged for up to four generations. At each passage, tumor tissues were cryopreserved (90% FBS and 10% DMSO), snap-frozen in liquid nitrogen for future use, and fixed in neutral buffered formalin for histological examination.

RNA extraction and quantitative real-time PCR

Total RNA was extracted from fresh PDX samples with TRIzol reagent (Invitrogen) in accordance with the manufacturer’s protocol. RNA concentration was quantified using a Nanodrop (Thermo Scientific, Hemel Hempstead, UK) and diluted to 100 ng/ml in RNase-free water before using the TransScript II One-Step gDNA Removal and cDNA Synthesis SuperMix kit (TransGen Biotech) according to the introduction. The cDNA was kept at − 20 °C until used for qPCR. Quantitative real-time RT-PCR analyses were performed using the SYBR^® Green Realtime PCR Master Mix (TOYOBO, BIOTECH CO., LTD), the BIO-RAD CFX96TM Real-Time System and Bio-Rad CFX Manager 2.1 software (Bio-Rad Laboratories, Inc., USA). The primer sequences for human FGF3 were: forward, 5′-ATGCTTCGGAGCACTACAGC-3′, and reverse, 5′-CCGTTCACAGACACGTACCA-3′. The primer sequences for human FGF4 were: forward, 5′-CTATGGCTCGCCCTTCTTCA-3′, and reverse, 5′-CCATTCTTGCTCAGGGCGAT-3′. The primer sequences for human FGF19 were: forward, 5′-AGATCAAGGCAGTCGCTCTG-3′, and reverse, 5′-GAGTACTGAAGCAGCCCCTG-3′. The primer sequences for human GAPDH were: forward, 5′-TTTGGTATCGTGGAAGGACT-3′, and reverse, 5′-AGTAGAGGCAGGGATGATGT-3′. The reaction conditions were an initial step of 95 °C for 60 s, followed by 40 cycles of denaturation at 95 °C for 15 s and annealing for 15 s, and extension at 72 °C for 45 s. Expression was normalized for RNA loading using GAPDH primers designed to span an intron and relative expression in each sample was calculated. To compare the difference of mRNA expression, fold changes were calculated by dividing the samples into two groups: CNV ≥ 5 or CNV < 5. Results are expressed as means ± SEMs. The semi-quantitative RT-PCR analyses from one individual experiment were repeated three times with comparable results.

H&E staining

Histopathology of primary P0 tumors and xenografts were evaluated using H&E staining according to standard method with an H&E staining kit (C0105, Beyotime, China). Results were reviewed by two independent pathologists.

Chemosensitivity of PDX model

Patient-derived xenografts models with 3 or more passages were used to evaluate chemosensitivity. Tumor bearing mice (~ 150 mm³ tumors) were randomized into two groups: treatment (paclitaxel, 10 mg/kg; platinum, 5 mg/kg, ip) and controls (saline). Treatment regimens were consistent with patient treatment and doses were determined from the literature [22–24]. Treatment was given once weekly for 3 cycles. Tumor size was measured twice weekly as described and tumor volume (V) was calculated as follows: V = L*W²/2 (L, length, long diameter of tumor; W, width, short diameter of tumor). All animal procedures were approved by the ethics committee for animal experiments at Peking University Cancer Hospital. Tumor growth inhibition (TGI) = (1 − ΔT/ΔC) × 100% (ΔT = tumor volume change of the drug-treated group on the final day of the study, ΔC = tumor volume change of the control group on the final day of the study).

Targeted next-generation sequencing and data analysis

Genomic DNA was extracted from P4 xenografts of each PDX model using a QIAamp DNA Mini Kit (QIAGEN Ltd., Crawley, UK) according to the manufacturer’s instructions. Extracted DNA was evaluated using a Qubit fluorometer (Invitrogen, Carlsbad, CA) and 1% agarose gel electrophoresis. A custom 483 cancer-related gene panel was used (Additional file 1: Table S1) [25]. The capture-based library was generated from 500 ng of DNA from each sample using a KAPA Hyper Prep Kit according the manufacture’s instruction (Kapa Biosystems, Boston, MA, USA), followed by Agilent’s SureSelectXT Target Enrichment System (Agilent Technologies, Santa Clara, CA). Library quality was assessed using Agilent 2100 Bioanalyzed on-chip electrophoresis (Agilent Technologies, Inc), and the library was quantified with an Agilent QPCR NGS Library Quantification Kit (Agilent Technologies, Inc). The library was sequenced on an Illumina Hiseq 2000 system (Illumina, San Diego, CA).

Sequencing adapters and low quality reads were filtered and the final Q20 and Q30 of all samples were > 90 and > 85%, respectively. The BWA software with default parameters was used to align sequencing reads to the human reference hg19 genome, and Picard was used to mark duplications. The aligned reads achieved coverage of > 99% of the target region with a mapping rate of > 95%. The average sequencing depth was more than 1200 × per sample. Mutations were called using Samtools, Mutect, and Varscan software, and were annotated using Annovar software. CNV analysis was performed using Event-Wise Testing algorithm according to previous reports [26]. A neutral copy number for each exon for each gene was established using normal lymphocyte samples.

Statistical analysis

Data were assessed for clinicopathological characteristics and transplantation success using a Chi squared test. An unpaired two-tailed Student’s t test was used to analyze the latency period of xenografts. Tumor growth between groups was compared using repeated-measured analysis of variance. A p value < 0.05 was considered statistically significant. Statistical analysis was performed with SPSS software for Windows, version 21 (SPSS Inc., Chicago, IL, USA).

Establishment of PDX models and patient clinical characteristics

A total of 188 esophageal carcinoma samples were obtained by endoscopic biopsy, and 40 PDX models were established in NOD/SCID mice at P1 (Fig. 1a; transplantation success 21.3%) Along with serial passaging, the transplantation rates from P1 to P2, P2 to P3, P3 to P4 were 70% (28/40), 89.3% (25/28), 100% (25/25), respectively. After the fourth generation, the PDX models became stable without further change in model formation, which could be used in future study. The overall transplantation rate for PDX models in our study was 13.3% (25/188). The median latency period of xenografts was decreased along with serial passage with 75.12 ± 19.87 days, 56.96 ± 17.17 days, 45.88 ± 21.38 days, and 28.04 ± 16.23 days at P1, P2, P3, and P4 passage (Fig. 1b, p < 0.05). Due to tissue limitations for the first two passages, a bio-specimen bank of xenografts from passage three were established, including formalin fixed paraffin-embedded blocks, viably cryopreserved and snap frozen tissues. Cryopreserved xenografts could be re-grown in mice, providing a renewable tissue source for future usage.

Histology of primary tumors and xenografts

Pathological assessment with H&E staining and analysis revealed consistent morphology and histology in xenografts compared with that of the corresponding primary tumor tissues (Fig. 1c; Additional file 1: Table S3). For case EC2, the PDXs kept poor differentiation from P1 to P5, which was in good accordance with the primary tumor. For case EC5, there was some variation in differentiation, with moderate differentiation in patient tumor and PDX of P1 to P6, except for P4, which was well differentiation (Fig. 1c). We then compared the histology type and differentiation of primary tumor and the two consecutive generation of P3, P4, and found that all the cases keep the same histology, and most cases kept the differentiation consistent (64%, 16/25), with some small changes in different passage (Additional file 1: Table S3).

Sensitivity of PDX models to chemotherapy

Five PDX models were treated with the clinically used chemotherapy agents, paclitaxel and cisplatin. Our data demonstrated that all of the PDX models had comparable therapeutic responses with that of the corresponding patients. The detailed data of the five PDX models and corresponding patients was shown in Fig. 2.

Genomic characteristics of ESCC PDXs

Using next generation sequencing (NGS) to molecularly profile xenografts, we identified 139 mutated genes from a total of 483 genes. Somatic mutations were variable (3–22; Fig. 3a, top). The observed mutation spectrum could be indicative of specific mutagenesis mechanisms occurring in ESCC. Among the 139 single-nucleotide variants, nucleotide transversion was over-represented, particularly C:G > T:A changes (Fig. 3a, bottom and right), which is consistent with ESCC patients in a previous study [27, 28].

Apart from the 139 genes with nonsynonymous single nucleotide variants (SNVs), 79 genes including insertions and deletions (InDels) in protein-coding regions were identified in the 23 PDXs. Thus, 191 non-silent somatic mutated genes were identified, with the range of 6-27 genes for each PDX (Fig. 3b, top). Of which, 96 genes were mutated in 2 or more samples. The top 10 mutated genes (frequency ≥ 7) in the 23 PDX models, and the 10 frequently mutated genes from other reports [27–32] appear in Fig. 3b. KMT2D was the most frequently mutated gene, involving the 7 frameshift inDel, one nonframeshift inDel and 7 missense mutation followed by EP300, PIK3CA, ALK and ERBB2. Frequencies of somatic mutations for each gene appear in Fig. 3b (right).

To interrogate the copy number alteration, 311 genes with copy number alterations (CNAs) were identified (median = 3.96, range 2.8–39.18). According to previous reports, genes with CNAs ≥ 5 are considered amplified [33, 34]. Finally, 63 genes were identified. And, we found 14 genes that were amplified in two or more PDX models (Fig. 3c). FGF3 (13/23), FGF4 (13/23), FGF19 (12/23). CCND1 (12/23) were the most frequently amplified genes and all had an amplification peak spanning 11q13.3. Additional peaks involving important cancer genes such as SOX2, MYC, VEGFA and CDK6 were also found. To demonstrate the amplification, we detected the expression of FGF3, FGF4, FGF19 using quantitative real-time PCR as reported papers [35]. All the samples were classified into two groups (CNA ≥ 5 and CNA < 5) according to the CNA tested by NGS, and CNAs ≥ 5 were considered amplified as it was mentioned above. As was shown in Additional file 2: Figure S1, in agreement with the CNA analysis, the relative mRNA expression of FGF3 (Additional file 2: Figure S1a, p < 0.001), FGF4 (Additional file 2: Figure S1b, p < 0.001) and FGF19 (Additional file 2: Figure S1c, p < 0.05) was enhanced in amplification group (CNA ≥ 5). The results was also in accordance with the previous studies [36, 37], indicating that most of the highly amplified genes showed elevated expression, and the CNAs provided a useful reference for amplification analysis.

Abnormal genes with SNV, InDel or amplifications were analyzed with DAVID Bioinformatics Resources 6.7, and significantly altered pathways (p < 0.05) were enriched and presented (Fig. 4a). The ErbB signaling pathway was the most significantly altered one, followed by the MAPK signaling pathway. Other pathways and abnormal genes appear in Fig. 4b. Most targets and pathways have been reported as candidate targets in carcinomas [7, 38–42]. These results suggested that the PDX models we established may be used for potential target selection and mechanism research.