Complex Patterns of Genomic Heterogeneity Identified in 42 Tumor Samples and ctDNA of a Pulmonary Atypical Carcinoid Patient

Tumor evolution underlies many of the most pressing challenges facing clinical precision oncology today, and a better understanding of this process has the potential to improve clinical care. We present a unique example of tumor evolution in an uncommon tumor type not previously featured in such studies, a pulmonary atypical carcinoid [an intermediate grade neuroendocrine tumor (NET)].

NETs arise from hormone-producing cells of the neuroendocrine system located throughout the body, and are highly heterogeneous, both genetically and pathologically. While once considered rare tumors (1) and thought to be indolent in nature (2), we now recognize that NETs have an age-adjusted incidence of 6.2 cases per 100,000 in the country where this study was undertaken, New Zealand (3), similar to the incidence of ovarian and cervical cancers (4). Pulmonary atypical carcinoids are mitotically active well-differentiated NETs, with one in five presenting with distant metastatic disease at diagnosis (5). They feature relatively few DNA variants, with recurrent variants occurring in chromatin remodeling genes such as MEN1 and ARID1A [reported in 25% and 10% of studied atypical carcinoids, respectively (6)]. In general, NETs feature few small DNA variants and often instead have large-scale chromosomal changes; for example, around 25% of pancreatic NETs lose a suite of 10 chromosomes, and a further 40% feature the loss of chromosome (Chr) 11 (7). Few genomic studies have been completed on pulmonary atypical carcinoids; however, comparative genomic hybridization studies identified recurrent deletions in Chr 11q (harboring MEN1) in around 60% of cases (8–10).

Tumor development is an evolutionary process with similarities to Darwinian natural selection, where tumor cells may be under multiple simultaneous “selective pressures” (11), including immune attack (12) and drug treatment (13). There are multiple debated models for the generation, propagation, and selection of variants throughout tumor development, including linear, branching, punctuated, and neutral evolution (11, 14), each of which may occur in different tumors or at different timepoints within the same patient. Tumor evolution has been associated with key clinical challenges, including cancer metastasis (15), immune evasion (12, 16), and drug resistance (13, 17, 18).

Sequencing the cell-free DNA (cfDNA) in a patient's blood plasma to identify ctDNA variants derived from tumor cells has been postulated to better represent the total disease burden (including potential tumor heterogeneity) than single tumor biopsies or resections (19, 20). However, we do not yet fully understand the detectability of ctDNA shed from different anatomic sites around a patient's body—and early results suggest that some tumor sites may be more easily detected than others (21).

In this study, we completed multimodal genomic analysis on samples collected at autopsy from 41 metastatic tumor sites from a single patient with an uncommon pulmonary atypical carcinoid, alongside a blood plasma sample and clinical biopsies, to catalog driver genomic alterations and hypothesize their order of accumulation. We compared the genomic variants identified in the tumors with those detected in the patient's blood plasma to consider how well the ctDNA analysis represented the genomic tumor heterogeneity, and in turn, infer metastatic sites that clinical ctDNA assays may poorly detect.

Informed written consent was obtained and a rapid research autopsy was completed in 8 hours after the patient's death under New Zealand Health and Disability Ethics committee approval 13/NTA/69/AM08, and in accordance with the Declaration of Helsinki. For a complete description of our ethical, legal, and logistical considerations, and tissue collection and sample processing protocol, please refer to Blenkiron and colleagues, 2019 (22).

DNA and RNA were extracted simultaneously using the Norgen formalin-fixed, paraffin-embedded (FFPE) RNA/DNA Purification Plus Kit (catalog no. 54300) from macro-dissected tumor regions of FFPE slides identified by a pathologist. DNA was extracted from blood buffy coat using the Qiagen QIAamp DNA kit (catalog no. 51304). Cell-free ctDNA was extracted from 5 mL defrosted blood plasma using the Qiagen QIAamp Circulating Nucleic Acid Kit (catalog no. 55114).

The 42 FFPE tumor samples, including 41 metastatic atypical carcinoid samples and one leiomyoma (in both cases classified by expert pathology review), were analyzed using whole-exome sequencing (WES) and transcriptome mRNA sequencing (RNA-seq). The two normal tissues were analyzed as germline controls by WES. Four tumors and one normal sample, overlapping with the WES and RNA-seq samples, were analyzed using 10X Chromium whole-genome sequencing (WGS). Clinical lung Biopsy 1 [a fine needle aspirate (FNA)], clinical subcutaneous breast Biopsy 2, and a blood plasma sample collected at autopsy were analyzed using a custom targeted DNA sequencing panel. Biopsy 1 was also analyzed by low coverage WES. All samples are described in Supplementary Table S1 and are named by the first two letters of the organ site.

DNA WES was completed as a service by Macrogen Inc., South Korea using at least 200 ng DNA per sample. Agilent SureSelect Human All Exon V6 capture (catalog no. 5190-8863) was used for 150 bp paired-end sequencing to a target depth of 200× coverage and processed through our standard exome variant calling pipeline (ref. 7; sequencing metrics in Supplementary Table S2). All variants were reviewed using Integrative Genomics Viewer (IGV, RRID:SCR_011793; ref. 23). Germline variants were curated using American College of Medical Genetics and Genomics guidelines (24). ADTEx (ref. 25; RRID:SCR_012059) was used for ploidy and zygosity estimations, while MuSiCa (26) and Signal (27) were used for mutational signature estimation, according to the tool creators’ instructions. ScarHRD was used to search for enrichment of homologous repair deficiency according to the tool creators’ instructions (28). Tumor mutational burden in mutations/Mb was calculated using the number of coding single-nucleotide variants (SNV) and indels and was compared with data from Supplementary Table S2 of Lawrence and colleagues, 2013 (29). Gene set enrichment analysis was completed using GenesetDB (30) with a FDR cutoff of 0.05. IQ-TREE (31) was used for phylogram creation from a FASTA alignment file of all variant sites, using ascertainment bias correction to account for constant sites not included in the alignment. Intertumoral heterogeneity and putative evolutionary pathways for the tumors of this patient were inferred using the LICHeE (32) and REVOLVER (33) software.

Exon target capture RNA-seq was completed as a service by Macrogen on the same tumors analyzed by DNA WES. A total of 500 ng RNA was used in the Illumina TruSeq RNA exome kit (catalog no. 20020189). Sequencing libraries underwent 100 bp paired-end sequencing with an expected output of 30 million reads/sample. Raw sequencing data were processed through our standard expression analysis and variant calling pipelines (ref. 7; sequencing metrics in Supplementary Table S3). Fusion analysis was completed using STAR-Fusion.

10x Chromium linked-read WGS was generated using high molecular weight genomic DNA extracted from fresh-frozen samples using the manufacturer's protocol (34) and the 10x Chromium Genome Library Kit & Gel Bead Kit v2 (catalog no. PN-120258). Sequencing libraries underwent 150 bp paired-end sequencing at a targeted output of over 200 Gb per sample on the Illumina HiSeq X Ten platform and were processed using the 10x “Long Ranger’ pipeline (RRID:SCR_018925, sequencing metrics in Supplementary Table S4).

A custom Thermo Fisher Scientific AmpliSeq HD panel was designed covering the 151 shared and private variants detected in WES (Supplementary Table S5). Sequencing libraries were prepared from 20 ng genomic DNA from tumor biopsies and cfDNA from blood plasma and sequenced on the Ion S5 platform with Ion 540 chip and 200 bp Ion Chef protocol (sequencing metrics in Supplementary Table S6). Data were analyzed using Ion Reporter 5.18 “AmpliSeq HD for Liquid Biopsy – w2.5 – DNA – Single Sample” pipeline with a hotspot file covering the known variant sites, and default filtering parameters.

Low-coverage WES of Biopsy 1 was completed as a service by Grafton Clinical Genomics. A total of 200 ng DNA was used in the Agilent SureSelect XT2 library kit with clinical research exome v1 capture probes, sequenced on the Illumina NextSeq500 to achieve a target of 40 million paired-end reads and processed through our standard exome variant calling pipeline (ref. 7; sequencing metrics in Supplementary Table S7).

Access to the datasets generated in this study is controlled out of respect for the deceased patient, however, may be requested through by a data management committee, made available from the European Genome-phenome Archive (https://www.ebi.ac.uk/ega/home; accession number EGAS00001006530).

A 69-year-old woman who had never smoked was diagnosed in 2007 with an atypical carcinoid tumor at the base of her right lung. Surgical resection was attempted but abandoned because of significant pulmonary artery involvement. When she died 10 years later, the tumor had spread to 90 anatomic locations detected on imaging, including brain, eyes, thyroid, liver, pancreas, kidneys, vertebrae, cranium, and many subcutaneous sites (clinical course summarized in Fig. 1A). A rapid research autopsy was completed within 8 hours after her death. A total of 358 tissue samples were collected from the 90 tumors. All tumors were histologically similar; therefore, a representative subset of tumors to be analyzed genomically was selected to include a diverse range of widely separated anatomic sites. We selected at least one sample of all large metastatic foci, and selected multiple samples for especially large lesions (16 cm³) to allow investigations of intratumor heterogeneity. A total of 42 tumor samples (including one uterine leiomyoma) and two non-tumor samples (heart muscle and blood buffy coat) underwent multiple genomic analyses (Fig. 1B). In addition, we included two tumor biopsies taken as part of clinical care (primary lung tumor FNA at diagnosis and breast tumor biopsy) and one blood plasma sample taken as part of a research project ten months before the patient passed away.

WES was performed on 41 metastatic atypical carcinoid samples, one leiomyoma and two normal samples (Fig. 1C). The benign uterine leiomyoma (Ut1) was, as expected, genomically distinct from the metastatic atypical carcinoid samples, so it will not be discussed further. Somatic chromosomal alterations were inferred from WES using ADTEx, revealing that all metastases sampled at autopsy shared the gain of one copy of Chr 5, loss of one copy of Chr 6, a complex “shattering” event of one copy of Chr 11 and loss of one copy of Chr 21 (Fig. 2). Biopsy 1 (diagnostic primary lung tumor FNA sampled 10 years earlier) featured the gain of Chr 5 and loss of Chr 21; however, there was insufficient evidence to definitively confirm the presence of changes to Chr 6 and Chr 11 in this historical biopsy (Supplementary Fig. S1) due to limited DNA mass and data quality (10.8 million reads were generated, resulting in only 15% of the exome footprint being covered with >30× coverage). The loss of one copy of a 246 kb region of Chr 12p was found in 15 of 41 tumors sampled at autopsy (Fig. 2D), disrupting ETV6 and resulting in loss of heterozygosity (LOH) of 360 genes, including tumor suppressor CDKN1B.

A feature of all tumors sampled at autopsy was the “shattered” Chr 11, indicative of chromothripsis. The precise junctions of 12 breakpoints on Chr 11 were not unequivocally identified by short-read sequence data but were revealed by linked-read sequencing generated on four representative tumor samples in comparison with one normal sample. Seven regions were lost from one allele of Chr 11, and the six remaining regions were reassembled in a new order, while the other allele of Chr 11 was unaffected (Fig. 3, with evidence of each breakpoint in Supplementary Fig. S2). Genes PLA2G16, RIC3, and C11orf87 were disrupted, and two possible fusion gene products were created (DLG2-CHKA and C11orf74-GRIA4; Fig. 3), while all other breakpoints fell in intergenic regions. The DLG2-CHKA fusion involved joining the first three (of 42) exons of DLGA to exons 8–14 of CHKA. It was absent from STAR-Fusion RNA-seq fusion analysis completed on the same tumor samples, suggesting the transcript may be reduced in abundance by nonsense-mediated decay. The C11orf74-GRIA4 fusion joined the last 12 exons of GRIA4 to the end of C11orf74 after seven (of eight) exons and was transcribed and detected in the RNA-seq in most tumors, with premature truncation shortly after the breakpoint. The lost regions of Chr 11 included the locus containing KMT2A, a frequent site of LOH in NETs (35), and gene set enrichment analysis of the genes lost revealed significant enrichment for the olfactory signaling, matrix metalloproteinases, and folate metabolism gene sets (Supplementary Table S8).

No pathogenic or likely pathogenic germline variants were detected, consistent with the patient's diagnosis in her seventies with no known family history of cancer. Tumor mutational burden (TMB) was around 1 variant per coding MB, in the lower third of TMB when compared with other well-studied cancer types (ref. 29; Fig. 4A). Our patient's tumors featured relatively consistent somatic TMB, ranging from 0.7 mutations per MB for some primary lung tumor samples (which may be representative of an earlier genomic state) to around 1 for Pa1, In1, and Lu9 (pancreas, abdominal wall nodule, and lung; Fig. 4B). The mutational signature making the greatest contribution was Signature 3, associated with BRCA1/BRCA2 inactivation (Supplementary Fig. S3A); however, it was only significantly enriched in two tumor samples [abdominal wall nodule In2 (60%) and pancreatic metastasis Pa2 (41%)] and did not reach the threshold for statistical significance in other samples. There were no BRCA1/BRCA2 or other homologous recombination DNA repair–related gene mutations identified in the WES, despite adequate depth of coverage. No other mutational signatures were significantly enriched. To robustly exclude homologous recombination DNA repair deficiency in this patient's metastatic tumors, we analyzed exome data from the tumors using the scarHRD package (28). All tumors, including the two that had a significantly enriched mutational signature 3 contribution (In2 and Pa2), had an HRD-sum score < 20; whereas known HR-deficient tumors invariably have HRD-sum scores of > 40 (ref. 28; Supplementary Fig. S3B).

All tumors sampled at autopsy shared 39 small variants (Fig. 5A). None of these shared variants were found in ClinVar, COSMIC or The Cancer Genome Atlas, and only a heterozygous variant in ARID1A was predicted to be a driver by the Cancer Genome Interpreter (CGI) algorithm (ref. 36; Supplementary Table S9). This single base pair deletion [chr1:g.27100383del, p.(Q1365Hfs*116)] in the “hypermethylated in cancer” domain of exon 17 causes a frameshift at amino acid residue 1365, leading to premature protein truncation.

A custom Thermo Fisher AmpliSeq HD sequencing panel covering all variants identified in tumors sampled at autopsy was used to interrogate the two clinical biopsy samples. Of the 39 variants shared by all tumors sampled at autopsy, just nine were detected in Biopsy 1 (diagnostic biopsy of lung primary tumor, 2007), sampled before there was evidence of metastatic disease by imaging or clinical examination, and are assumed to have accumulated earliest in tumor development (Fig. 5B). Of these, only three variants were predicted to be protein altering (in RRP12, MMP9, and MASTL; Supplementary Table S10); however, there is no evidence that these variants are associated with tumor development (36, 37). Those not detected in Biopsy 1 were sufficiently covered in the sequencing data to call as unequivocally absent in that biopsy (depth ranged between 1,065 and 9,463 unique DNA molecules). Surprisingly, the ARID1A putative driver variant was not detected despite the variant site being covered by 7,240 unique DNA molecules, so it likely developed after the initial diagnosis (Supplementary Fig. S4). Therefore, the gain of one copy of Chr 5 and the loss of one copy of Chr 21 (identified from low-coverage WES of Biopsy 1) preceded the development of the ARID1A variant.

A subset of samples from within the lung tumor carried very few additional variants over the 39 shared variants found in all lesions at autopsy (e.g., Lu5 and Lu7; Fig. 5B) and may possibly represent the most genomically similar tumors to the patient's historical primary lung tumor. Two variants were identified in all tumor samples situated outside the lung and three lung samples (GPAM and EPS8L2; Fig. 5A and B). Visual inspection of the sequencing data in IGV confirmed the absence of these variants in the remaining nine lung samples, despite adequate sequencing depth. Two further genomic groups were identified in the metastatic tumors sampled at autopsy. “Metastatic group 1,” identified by nine variants, spread to multiple vertebrae, the pancreas, thyroid, both kidneys, cranium, liver, subcutaneous sites, eye, and was also found in one sample from the lung tumor (Lu1, shown in group 2 in Fig. 5B due to subclonality). “Metastatic group 2,” identified by 11 variants, was detected in subcutaneous sites, pancreas, brain, vertebra and was also found in the lung tumor. Because all major genomic groups were represented within the lung tumors (Lu1 in “metastatic group 1,” Lu8 and Lu9 in “metastatic group 2,” apparent in Fig. 5C) it is possible that the metastatic clones may have arisen from this location. Within-sample genomic heterogeneity was evident in lung and pancreatic samples Pa1, Lu1, Lu8 and to a lesser extent Pa5, with evidence of both metastatic clones intermixed in these samples.

Biopsy 2 (subcutaneous breast, 2014) carried the same variants as sample Sc7 (from the same anatomic site sampled at autopsy, 2017), indicating that this lesion did not develop further genomic alterations in the 3 years before the patient died and suggesting that it remained genomically stable, despite its large size (19 cm³). Two rounds of chemotherapy were given in the interval between these two samples being taken (Fig. 1A); however, no additional variants consistent with selection driven by chemotherapy were identified.

Apart from the two clinical biopsies, all other samples were collected at a single timepoint (at autopsy). Nevertheless, we were able to apply two bioinformatic methods to this DNA dataset to first explore putative tumor phylogeny and subclonal structure (LICHeE) and secondly to infer the possible trajectory of variant acquisition (REVOLVER). The putative subclonal structure and evolutionary trajectory that appear most likely given this patient's metastatic genomic landscape are shown in Supplementary Fig. S5.

Having cataloged the shared and private small somatic variants held by this patient's sequenced tumors, we identified those detectable in a blood plasma cfDNA sample taken 10 months before the patient died. A custom Thermo Fisher Scientific AmpliSeq HD panel covering all the shared and private variants identified in tumors sampled at autopsy (Supplementary Table S5) was used for ctDNA sequencing, achieving a median depth of coverage greater than 6,000 unique molecules across variant sites.

The variants detected at the highest frequency in the ctDNA were shared by all tumors and, by inference, present in all cancer cells in the patient's body (at 20%–40% variant molecules; Fig. 6A). Two variants identifying metastatic group 2a (SLIT1 and TCF12) were detected at a similarly high level (around 20% variant molecules). Variants from metastatic groups 1 and 2b were detected at a markedly lower level (at 1%–3% variant molecules and 0.2%–0.6% variant molecules, respectively).

We searched for the presence of private variants, unique to each sampled site, in ctDNA to infer whether that lesion was uniquely detectable in the blood plasma, given that private variants would be shed from a single tumor site. (Fig. 6B). The only private variants detected were derived from subcutaneous breast sample Sc7 (four variants, each detected with 5–8 unique molecules), suggesting this may be a clear site of peripherally detected ctDNA release. No other private variants were detected with more than three unique molecules [the manufacturer's recommended cutoff for reporting variants (38), explored in Supplementary Fig. S6].

This study investigated a large number of metastatic tumors (two clinical biopsies, 42 tumor samples collected at autopsy and a blood sample collected during patient care) in a single patient with pulmonary atypical carcinoid cancer. Our aim was to provide a high-resolution understanding of the heterogeneity and putative developmental trajectory of this uncommon cancer type, complementing studies of more common cancer types that often involve large numbers of patients each with relatively small numbers of metastatic tumors. While a range of small variants and copy-number alterations have been reported previously in pulmonary atypical carcinoid tumors (9, 10, 39–41), no studies have explored the heterogeneity and relative timing of accumulation of multiple types of genomic alterations in this tumor's development.

No pathogenic germline variants were detected, and while the diagnostic biopsy contained numerous small variants not found in the patient's germline, none of these were predicted to be oncogenic drivers. However, large-scale chromosomal changes were present in the diagnostic biopsy and all other tumors of this patient (gain of one copy of Chr 5, loss of one copy of Chr 21). Given the clear absence of other driver mutations, we hypothesize that these chromosomal changes are likely to have driven, or at least strongly contributed to, initial tumor development. This chromosome 5 amplification is consistent with previous reports in atypical carcinoid tumors of recurrent gene-level copy-number changes on the p arm of Chr 5 (40).

In addition, all metastatic tumors collected at autopsy had loss of one copy of Chr 6 and chromothripsis of one copy of Chr 11. While the precise Chr 11 chromothripsis event described here is not reported in other studies, one atypical carcinoid cancer has been described with chromothripsis of Chr 11 and Chr 20 (42), and Chr 11 chromothripsis has been reported in other types of NET (43). In addition, deletions on the q arm of Chr 11 have been reported in between 50% and 72% of atypical carcinoid tumors (8–10), consistent with the general disruption of genes on chromosome 11 in this tumor type.

All metastatic tumors collected at autopsy also contained numerous small somatic variants not found in the diagnostic biopsy. However, only one of these small variants was predicted to drive tumor development; a heterozygous ARID1A variant. Somatic ARID1A variants, often heterozygous (44), have previously been found in atypical carcinoid tumors (6, 39) and have also been implicated in the pathogenesis of numerous other cancer types due to the role of ARID1A in the SWI/SNF chromatin remodeling complex (45), where ARID1A variants cause increased cell proliferation, migration, and invasion through a variety of mechanisms (46). The absence of the ARID1A variant in the diagnostic lung biopsy suggests that it most likely drove later tumor progression and was not an initiating driver variant, despite being identified as a driver in other cancer types and being recurrently detected in previous lung carcinoid studies (6, 39).

We can hypothesize a putative genomic evolution process (Fig. 7) that is consistent with the genomic heterogeneity observed at autopsy and in the two clinical biopsies. The early lung tumor represented by the diagnostic biopsy featured gain of Chr 5, loss of Chr 21 as well as nine small variants, none of which are known or predicted cancer drivers. Genomic features of the tumors collected at autopsy are consistent with further mutation after the lung biopsy was undertaken, followed by an evolutionary event such as a selective sweep, that fixed into the genomes of subsequent tumors a loss of function variant in ARID1A as well as 29 small non-driver variants and additional chromosomal alterations including Chr 11 chromothripsis and Chr 6 loss. There are plausible selective advantages to the tumor of each, including Chr 6 LOH of the HLA genes responsible for presenting peptides to T cells, hypothetically reducing the immune system's ability to recognize neoantigens presented by the tumor cells (47). The data are consistent with subsequent branching events, including the accumulation of a second somatic “hit” to EPS8L2 (where the first was somatic LOH on Chr 11), before two parallel lineages developed, metastatic group 1 and 2 (Fig. 7). In metastatic group 1, some variants of unknown significance deserve further study. For instance, truncation of tumor suppressor ETV6 [a common fusion gene partner in breast and thyroid cancer (48)], following the LOH of a region of Chr 12 may be functionally significant given this gene's role in development (49). The missense SLIT1 variant was the only additional protein-affecting variant carried by tumors in metastatic group 2 and may play a role in angiogenesis and migration (50). Interestingly, the two metastatic groups in the pancreas appear to have arisen by two independent seeding events. However, in the literature there is not overwhelming evidence for consistent genetic drivers of metastatic progression (51, 52) and we have no evidence in our data to suggest the additional genomic variants accumulated by each of the two metastatic groups (e.g., in EPS8L2, ETV6, and SLIT1) were necessary to drive their metastasis. Given a general propensity of cancers to metastasize to lung (53), it is possible that the presence of multiple metastatic groups in the lung may be the result of tertiary metastatic seeding events returning tumor cells that have evolved elsewhere to the site of the lung primary tumor.

The most likely subclonal structure and mutational progression hypothesized by the LICHeE and REVOLVER methods are consistent with the simple clustering analysis and phylograms shown in Fig. 5. While this consistency is interesting, we agree with the widely held opinion that inference of evolutionary history based on single timepoint data should be considered hypothetical (54).

ctDNA analysis is emerging as a valuable minimally invasive clinical test for patients with cancer (55, 56), capable of identifying the evolution of resistance mutations to targeted therapy (57). Here, we had an opportunity to compare the defined solid tumor DNA variant burden with that detectable circulating in the patient's blood plasma. The variants most readily detectable were shared by all or most tumor cells. The overall genomic heterogeneity was clearly detectable, with higher frequency variants being common to all or most lesions and lower frequency variants representing evolutionary outcrops. In the clinical setting where ctDNA sequencing would likely be completed without extensive solid tumor sequencing, ctDNA analysis would have identified the truncal variants that are shared by all metastases. While some other past studies have compared tumor variants identified in metastatic solid tumors and ctDNA (58, 59), few have investigated differences in the detection of ctDNA between anatomic sites. In this study, only the private variants in a single metastasis, breast lesion Sc7, was detected in this patient despite the larger size of other tumors including lung, thyroid, and pancreas (Fig. 6). The detection of ctDNA in the blood plasma derived from localized brain, thyroid, pancreas, and kidney tumors has been poor in some previous studies (21, 60, 61), as in our study, suggesting fundamental barriers to ctDNA detection including the blood–brain barrier, mucinous features (21), tumor desmoplasia (62) and first-pass metabolism by the liver (21). However, a past study in melanoma has identified that subcutaneous disease was not well represented in the ctDNA (63), in contrast to our finding. This suggests that the use of ctDNA for early detection of small treatment-resistant subclones may be challenging until more sensitive methods are developed and evaluated, especially if these subclones evolve in tissue types known to export little ctDNA (21, 60–62).

We acknowledge that studies that analyze many tumors from a single patient, while revealing detailed and valuable insights about that patient's cancer, cannot necessarily be generalized to other patients. Furthermore, all samples except for two clinical biopsies were collected at a single timepoint (at autopsy), limiting the extent of computational tumor evolution inferences that could be made, and the small size of these biopsies limited the scope of genomic assays that could be applied. Unfortunately, tissue of the two clinical biopsies was significantly depleted after clinical standard-of-care testing. This precluded additional analysis such as expanded targeted sequencing, which could more specifically place these samples in the evolutionary cascade of this patient's tumors. Furthermore, while RNA-seq data from this patient's tumors were useful for confirmation of expressed mutations and for identification of fusion gene products, a complex batch effect attributable to tissue processing artefacts in the autopsy samples could not be corrected using statistical approaches (Supplementary Table S3). This precluded the use of the RNA expression data in techniques such as molecular pathway analysis and differential expression analysis. Functional studies investigating the role of variants of uncertain significance, including chromosomal aberrations, were outside the scope of this study but provide an important avenue of future research.

In summary, DNA sequencing data generated from 42 widespread tumor samples, alongside clinical biopsies and blood plasma, have provided a high-resolution understanding of the genomic heterogeneity and putative evolution of this uncommon cancer type. In particular, this study provides evidence of the critical early role that chromosomal alterations may play in atypical carcinoid tumor progression. ctDNA analysis readily identified truncal variants but was unable to detect the private variants of most tumor subclones, despite high depth of coverage. This has implications for the use of ctDNA in cancer evolution studies and encourages further evaluation in expanded patient cohorts. We are grateful to the patient and her family for their contribution to scientific research.

No disclosures were reported.

T.J. Robb: Conceptualization, data curation, formal analysis, investigation, visualization, writing-original draft, writing-review and editing. P. Tsai: Formal analysis, validation. S. Fitzgerald: Investigation, writing-review and editing. P. Shields: Investigation, writing-review and editing. P.S. Houseman: Investigation. R. Patel: Investigation. V. Fan: Formal analysis. B. Curran: Formal analysis. R. Tse: Investigation, autopsy collections and pathological sample review. J. Ting: Investigation, pathological sample review. N. Kramer: Investigation, writing-review and editing, autopsy collections and protocol design, pathological sample review. B.J. Woodhouse: Investigation, project administration, writing-review and editing. E. Coats: Investigation, tissue sample cataloguing post-autopsy. P. Le Quesne Stabej: Data curation, writing-review and editing. J. Reeve: Investigation, radiological data curation. K. Parker: Investigation, project administration, writing-review and editing, patient and whānau communication/engagement. B. Lawrence: Conceptualization, funding acquisition, writing-review and editing, patient and whānau communication/engagement. C. Blenkiron: Conceptualization, supervision, methodology, writing-review and editing. C.G. Print: Conceptualization, supervision, funding acquisition, methodology, writing-review and editing.

T.J. Robb was supported by the Auckland Medical Research Foundation (1217001), and the wider project was funded by the Translational Medicine Trust. Figures 1C, 5C, and 6B were created with BioRender.com. The authors thank Grafton Clinical Genomics, led by Jason Copedo, and Auckland Genomics, led by Liam Williams, for their assistance with genomic assays, Nicola Dryden of Thermo Fisher Scientific for her assistance with optimizing the custom targeted sequencing panel, and the Auckland Regional Biobank for assistance processing blood plasma samples. Above all, the authors thank the patient and their family for their generous donation to scientific research.