A Diffusion-like Process Accommodates New Crypts During Clonal Expansion in Human Colonic Epithelium.

BACKGROUND & AIMS: Colorectal cancer (CRC) is thought to arise when the cumulative mutational burden within colonic crypts exceeds a certain threshold that leads to clonal expansion and ultimately neoplastic transformation. Therefore, quantification of the fixation and subsequent expansion of somatic mutations in normal epithelium is key to understanding colorectal cancer initiation. The aim of the present study was to determine how advantaged expansions can be accommodated in the human colon.
METHODS: Immunohistochemistry was used to visualize loss of the cancer driver KDM6A in formalin-fixed paraffin-embedded (FFPE) normal human colonic epithelium. Combining microscopy with neural network-based image analysis, we determined the frequencies of KDM6A-mutant crypts and fission/fusion intermediates as well as the spatial distribution of clones. Mathematical modeling then defined the dynamics of their fixation and expansion.
RESULTS: Interpretation of the age-related behavior of KDM6A-negative clones revealed significant competitive advantage in intracrypt dynamics as well as a 5-fold increase in crypt fission rate. This was not accompanied by an increase in crypt fusion. Mathematical modeling of crypt spacing identifies evidence for a crypt diffusion process. We define the threshold fission rate at which diffusion fails to accommodate new crypts, which can be exceeded by KRAS activating mutations.
CONCLUSIONS: Advantaged gene mutations in KDM6A expand dramatically by crypt fission but not fusion. The crypt diffusion process enables accommodation of the additional crypts up to a threshold value, beyond which polyp growth may occur. The fission rate associated with KRAS mutations offers a potential explanation for KRAS-initiated polyps.

Background and Context

As colorectal cancer is thought to arise from outgrowth of crypts harboring excess genomic alterations, quantification of mutation accumulation and spread within this tissue is key to understanding disease initiation.

New Findings

Exploiting loss of the cancer driver KDM6A, we demonstrate that newly generated crypts resulting from increased crypt fission are accommodated by mass movement of surrounding crypts in a diffusion-type process.

Limitations

Only 2 cancer driver genes, KDM6A and KRAS, are modeled here, which does not preclude the existence of alternative expansion mechanisms.

Impact

The threshold fission rate beyond which diffusion cannot accommodate newly generated crypts is calculable defining when identifiable pathologies may form.

It is widely recognized that many renewing epithelia acquire a substantial burden of cancer driver mutations while remaining apparently normal.1, 2, 3 In the human colon, development of neoplastic disease is thought to be driven by elevated rates of gland replication or fission. Most notably, loss of the tumor suppressor gene APC generates adenomas in this way.4, 5, 6, 7 Yet throughout life, normal crypts also undergo crypt replication at a low rate, which can be elevated by advantaged mutations.8, 9, 10, 11, 12, 13, 14, 15, 16, 17

There does not appear to be an increase in the net density of crypts, or of colonic epithelial area, with age. This raises the question: how are local clonal expansions arising from elevated fission rates accommodated? One explanation might lie in crypt fusion. This process has recently been described and could counteract the consequences of fission., However, it remains unclear if fusion is an independent stochastic process or if it is locally coregulated with fission. The latter possibility may be particularly relevant to pro-oncogenic mutations as fusions at the edge of mutant patches could enable effective local invasion of wild-type crypts with a high probability of subsequent displacement of wild-type cells.

KDM6A (UTX) is an X-linked gene encoding a histone demethylase that specifically targets di- and tri-methyl groups on lysine 27 of histone H3. Inactivating mutations and deletions of KDM6A have been identified in a variety of human cancers, including colon, bladder, prostate, and esophageal cancer.21, 22, 23, 24 KDM6A featured among the 127 significantly mutated genes in The Cancer Genome Atlas study that analyzed 3281 tumors derived from 12 cancer types. KDM6A mutations are infrequent in CRCs (<4% of all tumors, COSMIC database).

Here, in seeking additional cancer driver events that can be visualized as somatic clones, we identify loss of KDM6A as possessing advantage in both intracrypt fixation and subsequent expansion. The large multicrypt clones resulting from elevated rates of crypt fission are investigated to study the impact of expansion on crypt packing and the role of crypt fusion in relieving overcrowding. The increased fission rate within KDM6A− clones is not accompanied by an increase in crypt fusion, suggesting the 2 processes are driven by independent mechanisms and that fusion does not act to relieve local overcrowding. Instead, it is demonstrated that new crypts generated by fission can be accommodated by localized crypt diffusion up to a threshold beyond which hyperplastic and neoplastic lesions may form.

Materials and Methods

Human Tissue

Normal colon tissue samples were obtained from Addenbrooke’s Hospital Cambridge and Norfolk and Norwich University Hospital under full local research ethical committee approval (Documents 15/WA/0131 and 17/EE/0265, and 06/Q0108/307 and 08/H0304/85, respectively) according to UK Home Office regulations. The study included 273 individuals aged 13 to 93 years. Colectomy specimens were fixed in 10% neutral buffered formalin. From areas without macroscopically visible disease, mucosal sheets were removed and embedded en face in paraffin blocks. Sections were cut at 5-μm thickness onto charged glass slides.

Histochemistry

mPAS staining

This was performed as previously published.

Immunohistochemistry

Antibodies are listed in Supplementary Table 1. For standard sections, immunohistochemistry was performed as previously published. For laser capture microdissection, tissue was cut at 10-μm thickness onto UV-irradiated PEN membrane slides (Zeiss, Oberkochen, Germany). Heat-induced epitope retrieval was performed in citrate buffer in a water bath at 76 °C for 16 hours. Counterstaining with Mayer’s Hematoxylin was performed manually for 15 seconds followed by blueing in tap water for 1 minute.

Supplementary Table 1

Primary Antibodies Used for Immunohistochemistry

Antigen	Antibody	Supplier	Titer
KDM6A	HPA002111	HPA	1:100
KDM6A	#33510	CST	1:200
MAOA	SC-271123	Santa Cruz	1:200
STAG2	LS-B11284	LSBio	1:1000

Experimental Pathology

Data acquisition

Clones were scored by brightfield microscopy, followed by scanning of sections using a Leica (Wetzlar, Germany) Aperio AT2 scanner. The DeCryptICS image analysis tool developed by Edward Morrissey and Doran Khamis (https://github.com/MorrisseyLab/DeCryptICS) was used to count total number of crypts as well as Fusion or Fission (FUFI) forms per section. This was followed by manual quality control using QuPath including classification of FUFIs into mutant/mutant (M/M), mutant/wild-type (M/W), or entirely wild-type.

Quality control

Within the dataset, 2 individuals aged 37 years with extreme average patch sizes were identified as outliers with respect to that measure and not included in subsequent analyses of patch sizes, fusion rates, and newly generated crypts.

DNA Extraction From FFPE Tissue

Laser capture microdissection

Crypts were harvested into lids of 0.2-mm radius polymerase chain reaction (PCR) tubes using a Leica LMD7000 Laser Microdissection System; 10 μL of Proteinase K solution from the Arcturus PicoPure DNA Extraction Kit (ThermoFisher, Waltham, MA) was added followed by lysis at 65 °C for 3 hours and inactivation at 95 °C for 10 minutes.

Extraction from sections for KRAS sequencing

The QIAmp DNA FFPE Tissue Kit (Qiagen, Hilden, Germany) was used according to the supplier’s protocol.

Library Preparation and Sequencing

Primers, PCR reaction components, cycling conditions, and processing for amplification are described in Supplementary Table 2, Supplementary Table 3, Supplementary Table 4, Supplementary Table 5, Supplementary Table 6. Samples were barcoded using the Fast Start High Fidelity PCR System (Roche, Basel, Switzerland) according to the supplier’s protocol. After pooling and purification by Clean & Concentrator Kit (Zymo Research, Irvine, CA) and size selection by PippinBlue (Sage Science, Beverly, MA), samples were sequenced using 150-base pair (bp) paired-end sequencing with 10% PhiX in-house on the Illumina (San Diego, CA) MiSeq platform.

Supplementary Table 2

Primers Used for Amplification of FFPE to Assess Amplifiability

Primer	Sequence (5' – > 3')
MAOA_150_F1	ACCCATCAGTTACTCCTTCCC
MAOA_150_R1	GGGATTAAAGCTGGGAGTTTCT
MAOA_150_F2	TAGCAGGGCCTTGAATCTGT
MAOA_150_R2	GATAGTGCCCAGAGTCACCA
STAG2_150_F1	GGAGAAGAAGACACAGTTGGATG
STAG2_150_R1	TTCTGTGAGGCATTTAGGGAAAA
STAG2_150_F2	CCTATGCTCGCACAACTATGAG
STAG2_150_R2	GGAAGCCACACATCCTCTCT
CASD1_150_F1	ACCTGGAAACCCTATGCTCAA
CASD1_150_R1	TGCAGCTATACATGCCAACC
LOC_150_F1	TCGTCTGCTTCATCCTCCTC
LOC_150_R1	GCCTAACATGCTTGGACCAC
MAOA_250_F1	TGCAAGTCTTAGGTTGGTTGC
MAOA_250_R1	TCAGTAATGGGTCATGTGCAAA
MAOA_250_F2	AAGACATGTAGGGTTGGGGC
MAOA_250_R2	CAGAACACCCTGCTCTAACCT
STAG2_250_F1	GACTCTAAGGCCAGGTCAGG
STAG2_250_R1	GGAGGTGAGTTGTGGTGTCT
STAG2_250_F2	GCCTAATCATTCTCCCTGACCT
STAG2_250_R2	TGGTGTCAAAATCCATTCCCTC
CASD1_250_F1	GGTTAGAGGAAGACAAAAGTGGA
CASD1_250_R1	CCTCAGTCCACACTTTGATACAC
LOC_250_F1	AGCTTACCTCTTTGTCTCTTCCT
LOC_250_R1	CAACCTCAAAGTATCACGTGGA
MAOA_350_F1	TTCCTTCAGAAATTGAATCCTTG
MAOA_350_R1	CCTGGGAGAAAGCAAAATCA
MAOA_350_F2	TCCCGGAGTATCAGCAAAAG
MAOA_350_R2	CATGAGAGACCCCCAAACAC
STAG2_350_F1	TCCGAATATTTTTGGTGCATT
STAG2_350_R1	CAGAGCCTTGATGAGTGCTG
STAG2_350_F2	TCTGAAGGAATGCTATGGTATGAA
STAG2_350_R2	TTGTCAAGGGTCATAGACACAA
CASD1_350_F1	CTTTGGGAAGCTTTGCGTAAAA
CASD1_350_R1	CGATTCAGGAAGATGTAAGCCA
LOC_350_F1	TCAGGAATGATGGTCTACGTGA
LOC_350_R1	TCTCAGCTCTATTCCGTGAGT

Supplementary Table 3

Primers Used for Amplification of KDM6A

Primer	Sequence (5' – > 3')
1_F	ACACTGACGACATGGTTCTACACGCTTTCGGTGATGAGGAAA
1_R	TACGGTAGCAGAGACTTGGTCTCCGTACCTGTCCAGTCCG
2_F	ACACTGACGACATGGTTCTACATCTTTCAGGGCAATTAAAGCATT
2_R	TACGGTAGCAGAGACTTGGTCTACAACCTACCTTTAAACTAGACTCA
3_F	ACACTGACGACATGGTTCTACAGTACAATTGGACCATGGCCA
3_R	TACGGTAGCAGAGACTTGGTCTAGTGCAGAGGTATTACTACAACTT
4_F	ACACTGACGACATGGTTCTACACAGGATGCCATTAAATGCTACTT
4_R	TACGGTAGCAGAGACTTGGTCTTCTGGGGAAATATGTGGCTTT
5_F	ACACTGACGACATGGTTCTACAATGCTGTGTCACATCCTCCA
5_R	TACGGTAGCAGAGACTTGGTCTACTTGTTTGCTACCTCTACTCCT
6_F	ACACTGACGACATGGTTCTACATGACAGATGAGACCAACAGGA
6_R	TACGGTAGCAGAGACTTGGTCTCAGGCTGAGAGACGCTAGG
7_F	ACACTGACGACATGGTTCTACACTGCCTACAAACTCAGTCTCTG
7_R	TACGGTAGCAGAGACTTGGTCTCAGAAAAGGGTCCATTGGCC
8_F	ACACTGACGACATGGTTCTACATAACCGCACAAACCTGACCA
8_R	TACGGTAGCAGAGACTTGGTCTTCTCTCAAAGTGTATAAAACCCAGT
9_F	ACACTGACGACATGGTTCTACACGACCTCTCTCTTCCACTGG
9_R	TACGGTAGCAGAGACTTGGTCTAATGCCTTGTTGTCCACCTG
10_F	ACACTGACGACATGGTTCTACAGGCTGCTCTCAATCACCTCT
10_R	TACGGTAGCAGAGACTTGGTCTGCAGTGCTGTTAGGTGTCTC
11_F	ACACTGACGACATGGTTCTACAGAGACACCTAACAGCACTGC
11_R	TACGGTAGCAGAGACTTGGTCTTCCCATCAACAAGGCAGAGA
12_F	ACACTGACGACATGGTTCTACAGCCATTTCAACAGCAACACC
12_R	TACGGTAGCAGAGACTTGGTCTGGGGCTCTGAGATTCTTCCA
13_F	ACACTGACGACATGGTTCTACAGGAAGAATCTCAGAGCCCCA
13_R	TACGGTAGCAGAGACTTGGTCTCACACTAACCTGCATGCCTT
14_F	ACACTGACGACATGGTTCTACAATGGACTTGTGCAAATGCCTAGTAA
14_R	TACGGTAGCAGAGACTTGGTCTTGGAGGTGGACATTTATCCAACAA
15_F	ACACTGACGACATGGTTCTACATGTTTTCCTGAGATCTAACCACA
15_R	TACGGTAGCAGAGACTTGGTCTCAAGGCCACGTATTACTGTAACA
16_F	ACACTGACGACATGGTTCTACATGTAGAACACTAAACTAGACTGCT
16_R	TACGGTAGCAGAGACTTGGTCTACACAGTATTAGAAACATGCCTTTT
17_F	ACACTGACGACATGGTTCTACAGTTCTGGGAGGAGGAGGAAA
17_R	TACGGTAGCAGAGACTTGGTCTAGCACAGGATAACTCTTTGCA
18_F	ACACTGACGACATGGTTCTACAAAACTTCCACAGGTATTTGTAGC
18_R	TACGGTAGCAGAGACTTGGTCTCCAACATGGCTTAGAAGATTTCC
19_F	ACACTGACGACATGGTTCTACAGTGGAAGTTGCAGCTACATGA
19_R	TACGGTAGCAGAGACTTGGTCTTGCTCCCTGGAACTTTCATG
20_F	ACACTGACGACATGGTTCTACAACCGTGTGCTAACCAATTGC
20_R	TACGGTAGCAGAGACTTGGTCTACAAACCATTCACAGTCACCT
21_F	ACACTGACGACATGGTTCTACAGGAGCTTCTTAATGTAGTTGATCC
21_R	TACGGTAGCAGAGACTTGGTCTGCTGAATAAACCTATACACTGGAAC
22_F	ACACTGACGACATGGTTCTACACTAATGGGTTCTTGGTGGCC
22_R	TACGGTAGCAGAGACTTGGTCTTGAACCCAATGAACAGTGCC
23_F	ACACTGACGACATGGTTCTACAGCTGGTCACAAATAATTTCTCCC
23_R	TACGGTAGCAGAGACTTGGTCTTGAGCTGGTTCTTCTTTTGTCC
24_F	ACACTGACGACATGGTTCTACAACTTGGAAAACTTTGTGGTGCT
24_R	TACGGTAGCAGAGACTTGGTCTCACTGCTGCTTCATAACCCA

F, forward primer; Number, amplicon number; R, reverse primer.

NOTE. Sequence includes Fluidigm CS adapters.

Supplementary Table 4

Primers Used for Amplification of Area Around KRAS Codons 12 and 13 as Well as the Mimic Gene PITPNM2

Primer	Sequence (5' – > 3')
KRAS_1_F	ACACTGACGACATGGTTCTACATAAGGCCTGCTGAAAATGACT
KRAS_1_R	TACGGTAGCAGAGACTTGGTCTATGGTCCTGCACCAGTAATATG
MIMIC_1_F	ACACTGACGACATGGTTCTACATAGCACCCAGCCAGCTTG
MIMIC_1_R	TACGGTAGCAGAGACTTGGTCTATGACCACCCATGAAATATGAGCT
KRAS_2_F	ACACTGACGACATGGTTCTACAGGTGGAGTATTTGATAGTGTATTAACC
KRAS_2_R	TACGGTAGCAGAGACTTGGTCTTAGCTGTATCGTCAAGGCAC
MIMIC_2_F	ACACTGACGACATGGTTCTACACCTGCTGTTCCTACAAAGCTG
MIMIC_2_R	TACGGTAGCAGAGACTTGGTCTAGGCTTCTCCCGTCTAAGGA

Note. Mimic amplicons contained the sequence context of KRAS codons 12 and 13, including 5 preceding and 3 subsequent nucleotides (total synchronous sequence: GAGCT GGTGGC GTA) at the same position within the amplicon. Sequence includes Fluidigm CS adapters.

F, forward primer; Number, amplicon number; R, reverse primer.

Supplementary Table 5

Polymerase Chain Reaction Components

Gene	Primers	dNTPs	5X Phusion HF Reaction Buffer (New England BioLabs),	Phusion High-Fidelity DNA Polymerase	Water
KDM6A	1 μM	0.5 mM	5 μL	1 U	To 25 μL
KRAS	1.2 μM pair 1 0.8 μM pair 2	0.5 mM	5 μL	1 U	To 25 μL

Note. For KDM6A, pre-amplification was performed in multiplex groups (group 1: primers 1, 4, 7…group 2: 2, 5, 8… group 3: 3, 6, 9…). For KRAS, each sample was amplified in 2 independent duplex polymerase chain reactions containing a KRAS and corresponding mimic primer pair.

Supplementary Table 6

Polymerase Chain Reaction Cycling Conditions Used for Amplification of KDM6A and KRAS

Step	Temperature	Time	Cycles
Initial denaturation	95 °C	2 min	1
Denaturation	95 °C	10 s	35
Annealing	60 °C	10 s
Elongation	72 °C	15 s
Final Elongation	72 °C	5 min	1

Note. Polymerase chain reaction was followed by ExoSAP-IT enzyme treatment (2 μL of enzyme for 5 μL of sample, ThermoFisher) at 37 °C for 15 minutes and 15-minute inactivation at 80 °C and 1:10 dilution in DNA Suspension Buffer (Teknova). For KDM6A, samples were further amplified using the Fluidigm Access Array according to the supplier’s protocol.

Sequencing Data Analysis

Scripts are available under the following: https://github.com/kemp05/.

KDM6A

Amplicons were extracted by starting, finishing, and containing the expected sequence in the middle. Then, at every nucleotide position excluding the primer, the number of reads corresponding to the reference genome as well as those containing a base change were recorded. This enabled calculation of the noise at every position (Supplementary Table 7). Candidate mutations were identified when the mutant allele frequency (MAF) was either >4 times the mean of the noise at that position or >3.29 times the standard deviation at that position (P ≤ .001). True mutations were called if present in all samples originating from the same patch in serial sections but absent in all wild-type samples from the same sections.

Supplementary Table 7

Summary of Mutations Identified in KDM6A− Patches

Gender	Sample	Genomic coordinate	Exon	Ref	Alt	Alt freq	Maximum alt freq	Change
Male	1A	X:45078391	intron 18	A	G	85.8	100	Intronic
Male	1B	X:45078391	intron 18	A	G	72.9	100	Intronic
Male	2A	X:45078391	intron 18	A	G	82.9	100	Intronic
Male	2B	X:45078391	intron 18	A	G	85.3	100	Intronic
Male	3A	X:45085892	24	C	A	51	100	S1154>STOP
Male	3B	X:45085892	24	C	A	0	100	N/A
female	4A	X: 45089773	25	G	A	15.8	50	W1193>STOP
female	4B	X: 45089773	25	G	A	10.4	50	W1193>STOP

Note. Frequencies derived from frequency of mutant reads in samples. Of note, the mutation found in intron 18 was present in 2 independent patches. Furthermore, the mutations S1154∗ and W1193∗ have previously been found in cancer samples (COSMIC database).

Alt, alternative nucleotide -> this is the mutation; Change, predicted consequence of the mutation; freq., frequency; Ref, nucleotide in reference genome.

KRAS

Corresponding forward and reverse reads were combined using PANDAseq 2.11 with default options. Amplicons were extracted by beginning and ending with the expected primer sequence and correct overall length (± 3 bp tolerance).

For codon 12 and 13 mutation calling, reads containing the sequences corresponding to wild-type as well as all possible mutations (Supplementary Table 8) were extracted. This yielded an MAF for every possible mutation in all 4 amplicons for every sample and revealed the noise. KRAS mutations were called if (1) >1000 reads were obtained for both KRAS amplicons and at least 1 mimic amplicon, and (2) the MAF in both KRAS amplicons was > 0.1% (corresponding to at least 10 mutant reads) but found at background levels in the mimic amplicons. These criteria correspond to ≥1.96 standard deviations or a P value of < .025 for the noisiest nucleotide position (G12D). The actual MAF for a particular mutation was calculated by subtracting the mean allele frequency (noise).

Supplementary Table 8

Sequence Context Used to Extract Mutant Reads From Sequencing Data

Context	Nucleotide change	Amino acid change
GAGCTGGTGGCGTA	WT	WT
GAGCTGATGGCGTA	G>A	G12D
GAGCTGCTGGCGTA	G>C	G12A
GAGCTGTTGGCGTA	G>T	G12V
GAGCTAGTGGCGTA	G>A	G12S
GAGCTCGTGGCGTA	G>C	G12R
GAGCTTGTGGCGTA	G>T	G12C
GAGCTGGTGACGTA	G>A	G13D
GAGCTGGTGCCGTA	G>C	G13A
GAGCTGGTGTCGTA	G>T	G13V
GAGCTGGTAGCGTA	G>A	G13S
GAGCTGGTCGCGTA	G>C	G13R
GAGCTGGTTGCGTA	G>T	G13C
GAGCTGGAGGCGTA	T>A	WT
GAGCTGGCGGCGTA	T>C	WT
GAGCTGGGGGCGTA	T>G	WT
GAGCTGGTGGAGTA	G>T	WT
GAGCTGGTGGGGTA	G>C	WT
GAGCTGGTGGTGTA	T>C	WT

Note. For each variant, the mutant allele frequency was calculated by dividing by the total read number.

Mathematical Modeling

Stem cell dynamics and crypt fission

The stem cell dynamics and fission rate associated with loss of KDM6A were mathematically modeled as previously described.

Crypt fusion and diffusion

Crypt fusion was modeled as a process parallel to fission, with the same duration. Therefore, the following applied:

All M/W FUFIs were considered fusions, but M/M FUFIs could be fissions or fusions. Therefore, the number of fission events and fusion events required for the preceding equation were not directly measurable. However, they were calculable by sampling FUFIs from the edge of mutant patches, which enabled calculation of the proportion of M/M FUFIs that are fusions (termed chi: χ) by accounting for the distribution of the mutational state of neighboring crypts (see supplemental materials for mathematical details).

To infer a diffusion coefficient, growth of a patch through initial mutation and subsequent fissions was modeled as a stochastically firing point source of mass at the clone centroid. Potential trajectories from mutation hit to patch of size 10 (i.e, stochastic event times) were simulated, and an ensemble diffusion coefficient was inferred by randomly drawing paths from the set of simulated potential trajectories. The inferred diffusion coefficient was then used in a theoretical study of patch expansion (see supplemental materials for mathematical details).

KRAS expansion

The data obtained here were combined with our previously published dataset and analyzed as described there.

Results

KDM6A-negative Clones Are Advantaged in Stem Cell Competition

We have previously used visualization of loss of X-linked genes as clonal marks to quantify human colonic stem cell dynamics. In attempting to expand this methodology to X-linked genes with cancer association, clonal loss of KDM6A was identified by immunohistochemistry with 2 independent antibodies on normal human colonic epithelia (Supplementary Figure 1A and B). Intracrypt dynamics that describe the accumulation of clones wholly populating entire crypts (WPC) from partly populated (PPC) transition forms were determined for KDM6A- clones using colonic FFPE sections from 120 patients aged 21 to 93 as previously described (Figure 1A and B). This revealed that loss of KDM6A confers a competitive advantage to affected stem cells shown as a decrease in the fraction of PPC supporting the accumulation of WPC (Figure 1C).

Supplementary Figure 1

Detection and laser capture of KDM6A− clones. (A) Representative images of immunohistochemistry (IHC) for cancer-associated X-linked genes on human colonic sections. Of the 7 assessed genes, a staining pattern indicative of loss-of-function mutation was only found for KDM6A (middle). (B) Serial human colonic sections stained with 2 independent antibodies against KDM6A: (i) human protein atlas, (ii) Cell Signaling Technology. KDM6A− crypts are circled and enlarged. (C) For assessment of optimal amplicon size, DNA extracted from FFPE sections was diluted to equivalents of 5, 20, and 100 crypts and used for PCR reactions (primers in Supplementary Table 2). Image of gel electrophoresis of PCR products on FFPE DNA diluted to equivalents of varying amounts of crypts. This revealed that an amplicon size of 150 bp ensures reliable amplification. (D) Image of IHC for KDM6A on laser capture slide (i) before and (ii) after laser capture microdissection. (E) Box plot showing reads obtained for each of the 24 amplicons used for sequencing of KDM6A− and control patches. Amplicon libraries from 4 patches failed quality control. Of the remaining 7 patches, mutations were identified in 4 at mutant allele frequencies consistent with patient sex and estimated stromal fraction within the captured material (Supplementary Table 4). These results both confirm antibody specificity for KDM6A and the clonality of multicrypt patches that share only a single mutation.

Figure 1

Detection of clonal loss of KDM6A in normal human colon. (A) Representative images of KDM6A− WPC (i) and PPC (ii). (B) Top: Frequency plot showing age-related behaviors of KDM6A− WPC (circles, darker) and PPC (squares, lighter). Red line: regression analysis showing accumulation of WPC (slope of accumulation of fixed clones [ΔC]: 6.04 × 10−6 per year) and 95% CI in gray. Bottom: PPC only on expanded y-axis. (C) Plot showing increased ratio of ΔCfix/frequency of partial clones (Cpart) for KDM6A (0.23, 95% CI 0.16–0.34) as compared with the neutral marks mPAS and MAOA (replotted from Nicholson et al). Error bars = 95% CI. (D) Representative image of large KDM6A− multicrypt patch highlighted by dashed line. (E) KDM6A cDNA structure annotated with sequenced areas (yellow), mutations found in COSMIC and mutations identified in KDM6A− patches (red).

KDM6A-negative Clones Expand by 5-Fold Increased Crypt Fission

Expansion of individual KDM6A− clones was recognizable as large patches that frequently exceeded 10 crypts (Figure 1D). To confirm the clonal origin of such patches, we used laser capture microdissection followed by targeted sequencing that covered 3.6 kb of exonic and flanking intronic sequence of KDM6A (24 amplicons), including sites frequently mutated in human cancers. No patch was found to carry more than 1 KDM6A mutation and mutant allele frequencies were in line with predictions from patient sex and stromal content, supporting clonality (Figure 1E, Supplementary Figure 1C–E).

The age-related size distribution of multicrypt clones was analyzed to infer the fission rate associated with KDM6A loss and compare it with those previously described for neutral (MAOA and mPAS) and advantaged clonal marks (STAG2). This revealed an age-related increase in the frequency of large clones (≥10 crypts/patch) for STAG2− and KDM6A− but not for neutral marks (Figure 2A). Loss of KDM6A generates a higher proportion of large patches than the other clonal marks, whereas STAG2 loss generates more clones due to a higher event rate (Figure 2B and C). From these patch size distributions, the crypt fission rate associated with loss of KDM6A was calculated to be to 3.6% per year (95% confidence interval [CI] 3.2–4.1), approximately 5-fold higher than the background homeostatic rate previously derived from neutral clonal marks (Figure 2D). Consequently, in individuals older than 80 years, 13.5% of KDM6A− clones are found as patches comprising more than 5 crypts compared with 4.8% for STAG2, 1.7% for mPAS, and 0.8% for MAOA.

Figure 2

KDM6A− patches lack significant overcrowding despite increased fission. (A) Plot showing mean frequency of large (≥10 crypts) patches for age groups shown. (B) Histogram showing the frequency of different patch sizes for mPAS, MAOA, KDM6A, and STAG2 across all ages. Inset shows patch size ≥2 crypts on expanded y-axis. (C) Dot plot of mean clone frequency plotted against mean average patch size for multicrypt clones in age groups shown for mPAS, MAOA, STAG2, and KDM6A. (D) Plot showing inferred fission rate/crypt/year for KDM6A compared with mPAS, MAOA, and STAG2 (data replotted from Nicholson et al). Error bars = 95% CI. (E) Simulation data showing the time in years taken for transitions between patch sizes 1–2 and 10–11 for (i) mPAS, (ii) STAG2, and (iii) KDM6A. Insets show 10–11 transition on expanded y-axis. (F) Image showing patch selection for assessment of crypt packing density within KDM6A− (blue) and adjacent wild-type crypts (brown). Inset indicates placement of borders between crypts. (G) Dot plots showing stromal fraction for STAG2− and KDM6A− and adjacent wild-type patches comprising 10 crypts.

KDM6A− and STAG2− Patches Lack Significant Local Overcrowding

Expanding clones are more likely to undergo additional fissions as the probability scales with the number of crypts present. Consequently, the interval between fissions decreases rapidly with increasing patch size. For example, mathematical modeling of the expansion of KDM6A− crypts indicates the median time taken to grow from 1 to 2 crypts is 19 years, but only 2 years to grow from 10 to 11 (Figure 2E). Therefore, recently formed larger patches might be expected to demonstrate overcrowding.

To test if larger patches are more densely packed, the area occupied by crypts and their surrounding stroma was determined for 24 STAG2− and 20 KDM6A− clones containing 10 crypts, the largest size for which sufficient data could be obtained (Figure 2F and G). The fraction of each patch occupied by stroma was then calculated. Adjacent to each mutant clone, 3 random groups of 10 crypts were defined as control “patches” and similarly analyzed (totaling 72 control patches for STAG2 and 60 control patches for KDM6A). Comparing the fraction of each patch occupied by stroma to adjacent wild-type groupings indicated a slight trend toward increased packing density for STAG2− as well as KDM6A− crypts that failed to reach significance (Figure 2G). Considering that a lack of overcrowding may stem from a decrease in crypt size, the areas of crypts were measured. This revealed that KDM6A− crypt sections are approximately 1.3 times the size of adjacent wild-type crypts (P < .001) (Supplementary Figure 2). No difference was found between STAG2− crypts and their wild-type neighbors (Supplementary Figure 2). Therefore, lack of overcrowding cannot be attributed to reduced crypt size for either STAG2 or KDM6A loss. We also considered the possibility of accommodation of clonal expansions by “squashing” of neighboring crypts. However, analysis of crypt eccentricities provided no evidence for the predicted flattening of crypt architecture that would result (Supplementary Figure 3, Supplementary Figure 4, Supplementary Figure 5).

Supplementary Figure 2

Size measurement of STAG2− and KDM6A− crypts. Crypt areas in multicrypt clones containing 10 crypts that are STAG2− or KDM6A− (24 and 20, respectively). The area of individual crypts was determined using QuPath software. As controls, adjacent or nearby wild-type groups of 10 crypts were analyzed in the same way. Jitter plots show areas of crypts. Black line = median, ∗∗∗∗P < .001.

Supplementary Figure 3

Eccentricities of crypts around STAG2− clones of size 10 as a function of radial distance from crypt centroid. Crypt squashing due to clone growth would impact wild-type (WT) (orange) crypts closest to the expanding patch, predicting a right to left decline in the crypt eccentricities in the plots above with distance from mutant (brown) crypts. This pattern is not observed.

Supplementary Figure 4

Eccentricities of crypts around KDM6A− clones of size 10 as a function of radial distance from crypt centroid. Crypt squashing due to clone growth would impact wild-type (WT) (light turquoise) crypts closest to the expanding patch, predicting a right to left decline in the crypt eccentricities in the plots above with distance from mutant (teal) crypts. This pattern is not observed.

Supplementary Figure 5

Comparison of mean eccentricities of mutant and adjacent WT crypts in clones. For each assessed patch of 10 KDM6A− or STAG2− crypts, the mean eccentricity of the mutant crypts is plotted on the x-axis (0 = perfect circle) against the mean eccentricity of the surrounding wild-type crypts on the y-axis. The diagonal line represents x = y, where the mean eccentricities of mutant and surrounding WT crypts correspond. Error bars = 95% CI around the mean.

Together these observations suggest that crypts even within relatively recent clonal expansions avoid overcrowding to largely achieve ambient density.

Evidence for Crypt Fusion

The lack of overcrowding in KDM6A− clones suggests a mechanism counteracting the localized increase in fission. An opposing process of crypt fusion has been recognized in mouse intestine. A homeostatic human fusion rate has been estimated by assuming equivalence in the rate of both fission and fusion. On a tissue-wide basis such a balance of rates could act to maintain constant crypt density. However, local advantaged expansions can be balanced only if fission and fusion are locally coordinated.

We first sought confirmation that fusion occurs. The evidence in human epithelium is based on identification of branched crypts within which clonal loss of mitochondrial CCO activity is restricted to one branch. These are interpreted as transition intermediates in an active fusion process. Analysis of en face tissue sections stained to visualize mPAS positivity confirmed the existence of rare heterotypic branched forms in normal human colonic epithelium. Analysis of more than 2 × 106 crypts in sections from 80 individuals containing mPAS+ clones identified 32 candidate mPAS+ branched forms that were either mixed (mutant and wild-type, M/W) or fully mutant (M/M) (Figure 3A). Of the 13 M/W forms, the positive epithelium was always restricted to 1 branch.

Figure 3

Crypt fission and fusion are independently regulated processes. (A) Schematic and representative mPAS-stained images of 3 types of fusion or fission forms. (B) Schematic representation of alternative origins of M/W forms. Hypothesis 1 (Hyp1): stem cell mutation in one branch of intermediate fission form followed by monoclonal conversion. Hypothesis 2 (Hyp 2): fission of a preexisting partially populated crypt with segregation of mutant and wild-type epithelium into each branch, followed by monoclonal conversion. (C) Comparison of M/W event frequencies simulated for hypotheses described in (B) and the observed frequency. Error bars = 95% CI. (D) Bar graphs showing numbers of FUFI types scored for (i) mPAS, (ii) STAG2, and (iii) KDM6A. (E) Representative images of FUFIs at patch borders: (i) KDM6A M/M FUFI, (ii) STAG2 M/W FUFI. Insets show enlarged FUFIs. (F) Schematic showing scoring of FUFI neighbors used for calculation of χ. (G) Dot plot showing neighboring crypt status of patch border FUFIs (M/M or M/W) for mPAS, STAG2, and KDM6A. Each dot represents 1 or more FUFIs for which neighbors were scored. (H) Plot comparing the value of χ for mPAS, STAG2, and KDM6A. (I) Plot comparing the derived fission and fusion rates for mPAS, STAG2, and KDM6A. mPAS and STAG2 fission rates correspond to data from Nicholson et al, replotted. Error bars = 95% CI.

An alternative interpretation is that branched crypts are intermediate fission forms whereby heterotypic staining arises due to mutations occurring or segregating into a single branch (Figure 3B). We formally considered this possibility using the fusion duration estimate derived by Baker and colleagues as well as our previous estimates of de novo mutation probability and clone fixation rates that together determine the frequency of monoclonal crypts present in individuals of different age. Within the relatively small number of branched crypts present, none are predicted to contain monoclonal crypt branches by either mechanism (Figure 3C). This suggests that heterotypic forms represent genuine intermediates in an active fusion process. Because the bulk of branched crypts are unstained and can represent intermediates in either fusion or fission we propose the agnostic term FUFI to describe these transition forms (Supplementary Figure 6).

Supplementary Figure 6

Defining the FUFI. (A) Representative image of KDM6A-positive bifurcating crypts sectioned longitudinally. (B) A FUFI is defined as 2 adjoined crypts viewed in a transverse section with 2 clearly discernible lumina lacking any separating gap.

Crypt Fission and Fusion Are Regulated Independently

To calculate crypt fusion rates, heterotypic and homotypic FUFIs were also evaluated for STAG2 and KDM6A loss by scoring approximately 3.9 × 106 and 1.8 × 106 crypts from 53 and 102 individuals, respectively. In total, more than 28,000 FUFIs were evaluated. This identified 151 and 63 FUFIs with STAG2 (18,928 clones analyzed) and KDM6A loss (5,353 clones analyzed), respectively (Figure 3D). These could be found as single events or within multicrypt clones (Figure 3E).

Assuming equal duration for fission and fusion (ie, the time window during which FUFIs are detectable), the fusion rate is accessible by proportionality. Specifically, the ratio of the frequency of observed fission FUFIs to the fission rate (independently calculated from the patch size distribution) would equal the ratio of the frequency of observed fusion FUFIs to the fusion rate. However, although all M/W FUFIs are considered fusions, M/M FUFIs can be either fissions or fusions. Therefore, to exploit this proportionality, the relative contribution of M/M FUFI events to the total number of fusions (termed chi: χ) needs to be determined. The value for χ is calculable on the basis that M/M FUFIs have a probability of being a fusion event that depends on the number of M and W neighbors present at the onset of the fusion process. Therefore, the status of crypts neighboring FUFIs at the patch border of multicrypt clones were scored as W or M (totaling 4, 121, and 49 for mPAS, STAG2, and KDM6A, respectively) (Figure 3F and G). Single FUFIs had W neighbors only. Averaging the M and W neighbors of M/M FUFIs revealed χ to be 0.03 for mPAS. As neutral marks such as mPAS generate mostly small clones, most fusion events are W/M leading to low values of χ. At the edge of larger expansions generated by advantaged marks, M/M fusions occur more readily. Correspondingly, the χ values for STAG2 and KDM6A are 0.26 and 0.35, respectively (Figure 3H). Rates of fusion using these values of χ were then estimated using the proportionality described previously (see supplemental information for details).

This analysis indicates similar crypt fusion rates for mPAS, STAG2 and KDM6A of 0.3% per year (95% CI 0.1–0.6), 0.4% (95% CI 0.3–0.7), and 0.7% (95% CI 0.3–1.4) (Figure 3I). Comparison of mPAS fission (0.7% per year; 95% CI 0.5–0.9) and fusion rates show that these closely correspond. This suggests that in homeostasis the rates of both processes are balanced and will act together to maintain constant crypt numbers across the tissue, as has been suggested previously. However, for mutations causing elevated fission rates there appears to be no evidence for a compensatory increase of the fusion rate. These analyses suggest that fission and fusion are independent processes and not coordinately regulated.

Crypt Diffusion Accommodates New Crypts Throughout Life

A striking feature of larger patches is that mutant crypts have over decades populated the territory initially occupied by multiple independent crypts without a significant increase in crypt density. In the absence of appreciable crypt fusion this suggests local adjustments to disperse crypts from the growing focus. With this rationale, we considered the possibility of random crypt movement in the form of a diffusion process.

In the colonic epithelium it is the crypts that are being diffused in the “space” of the surrounding stroma (Figure 4A). The diffusion coefficient (change in area per unit time) can be estimated based on crypt packing (measured in terms of crypt area per unit area of mucosa) and consideration of all possible sequences of fission events, constrained by the patient age and calculated using the mutation and fission rates.

Figure 4

Evidence for a crypt diffusion process. (A) Representation of the proposed crypt diffusion process. Fission generates higher local density which is relieved by diffusion. (B) Schematic showing crypt cross section with area of surrounding stroma defining a crypt domain. (C) Representation of areas measured to assess packing of mutant patch and surrounding crypts in rolling windows placed at different distances (r) from centroid of mutant patch. (D) Examples of radial variation in crypt packing in patches of 10 crypts when moving from mutant clone to 3 adjacent control groupings for (i) STAG2 and (ii) KDM6A. Points are data derived, the black line (gray ribbon) is the median (95% CI) theoretical stromal fraction as fitted from the diffusion model. Dashed line shows the average of the 25 most likely trajectories from initial mutation to clone of size 10, based on population average diffusion and neighborhood ambient stromal fraction. (E) As (D) but with rolling window data. (F) Density plot of values obtained for the diffusion coefficient in human colonic epithelium. (G) Representation of area of crypts affected if space is decreased by 1% or 5% surrounding a mutant patch of 10. (Crypt numbers not to scale.) (H) Boxplot showing simulated numbers of crypts affected if spacing is decreased by 1% or 5% respectively by addition of a patch of 10 mutant crypts. (I) Dot plot showing median frequency of newly generated crypts for mPAS, MAOA, STAG2, and KDM6A across 4 age bins. (J) Line graphs showing the stromal fraction resulting from crypt diffusion and crypt fission rates corresponding to the homeostatic (WT) rate as well as those associated with STAG2, KDM6A, and KRAS and multiples of 11 and 12 of the WT rate. Dotted line = whitespace fraction calculated from optimal hexagonal packing of circles. r = distance from centroid of patch in crypt domains. (K) Line graph showing results from simulations to calculate the cumulative probability of a clone developing into a lesion, defined as reaching the whitespace fraction for hexagonal packing shown in (J), over time. Plots using the fission rate associated with KRAS (17 times WT, none) are shown alongside those including therapeutic crypt fission inhibition applied for a decade either immediately after mutation acquisition (0–10), after 10 years (10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20) or 20 years (20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30).

To find evidence supporting a diffusion-type process we revisited the STAG2− and KDM6A− patches of 10 crypts and surrounding control patches. For each mutant clone (24 for STAG2 and 20 for KDM6A) and the corresponding 3 arbitrary control groups of 10 adjacent wild-type crypts, each crypt was spatially mapped in X/Y coordinates. Total patch areas were divided in crypt domains that were defined as the area occupied by a crypt and its share of surrounding stroma (Figure 4B). Areas of individual crypts as well as total patch area were measured and the distance between mutant clones and wild-type patch centroids, r, was determined (Figure 4C).

Using known mutation and fission rates, the potential trajectories from initial mutation to 10-crypt mutant clones were simulated and the most likely trajectories were used to calculate the age of the clone (Supplementary Figures 7 and 8). The diffusion coefficient was inferred to define an overall tissue process that best fits the differences in stromal densities between mutant and control patches. It describes how the burden of decreased stromal fraction resulting from a mutant clone is dispersed into the surrounding tissue over time. For older clones, the expectation is for the system to be close to the ambient density, whereas for younger clones, the perturbation to the local stromal fraction may still be evident. Examples of both presumptive young and old clones were readily detectable (Figure 4D, Supplementary Figures 9 and 10).

Supplementary Figure 7

Inferred ages of STAG2− clones. Theoretical distributions for the age of STAG2− clones in years covering the period from single crypt to a clone comprising 10 crypts. The theoretical density of patch age is calculated using the mutation rate, the monoclonal accumulation rate, and the fission rate (modeling fission as a stochastic birth process). Dashed lines show the most likely patch age for each clone derived from the average of the 25 most likely trajectories.

Supplementary Figure 8

Inferred ages of KDM6A− clones. Theoretical distributions for the age of KDM6A− clones in years covering the period from single crypt to a clone comprising 10 crypts. The theoretical density of patch age is calculated using the mutation rate, the monoclonal accumulation rate, and the fission rate (modeling fission as a stochastic birth process). Dashed lines show the most likely patch age for each clone derived from the average of the 25 most likely trajectories.

Supplementary Figure 9

Diffusion model prediction of stromal fraction changes for STAG2− clones.Plots show the stromal fraction measured for STAG2− as well as surrounding wild-type (WT) patches of 10 crypts. The radial distance r (in crypt domains) is measured from adjacent patch centroid to mutant patch centroid. The black line and gray ribbon represent the median and 95% CI theoretical stromal fraction as fitted from the diffusion model. Dashed line shows the solution from averaging the 25 most likely trajectories from initial mutation to clone of size 10, assuming population average diffusion and neighborhood ambient stromal fraction. Green box: patches for which a “rolling window” was applied. For each mutant patch, surrounding measurements included 2 areas comprising 3 mutant and 7 WT, 2 mutant and 8 WT, and 1 mutant and 9 WT crypts as well as 5 surrounding WT patches at varying distances (combined: 1 data point for the mutant patch and 11 for surrounding patches).

Supplementary Figure 10

Diffusion model prediction of stromal fraction changes for KDM6A− clones. Plots show the stromal fraction measured for KDM6A− as well as surrounding wild-type (WT) patches of 10 crypts. The radial distance r (in crypt domains) is measured from adjacent patch centroid to mutant patch centroid. The black line and gray ribbon represent the median and 95% CI theoretical stromal fraction as fitted from the diffusion model. Dashed line shows the solution from averaging the 25 most likely trajectories from initial mutation to clone of size 10, assuming population average diffusion and neighborhood ambient stromal fraction. Green box: patches for which a “rolling window” was applied. For each mutant patch, surrounding measurements included 2 areas comprising 3 mutant and 7 WT, 2 mutant and 8 WT, and 1 mutant and 9 WT crypts as well as 5 surrounding WT patches at varying distances (combined: 1 data point for the mutant patch and 11 for surrounding patches).

For a subset of 7 patches (4 for KDM6A, 3 for STAG2), a more detailed rolling window analysis was performed in which the preceding approach was applied but where fields of 10 crypts were moved outward from the mutant clone. Again, evidence of a reduction in stromal fraction consistent with perturbation in younger clones was observed (Figure 4E, Supplementary Figures 9 and 10). The diffusion coefficient was found to be 1.05 crypt domain areas per year (95% CI 0.339–9.70) (Figure 4F). Reassuringly, testing a null hypothesis, that there is no diffusion-type process and therefore no radial dependence in crypt packing, by considering neighborhood ambient densities of crypts across all patches (mutant, wild-type, or mixed), generated a significantly worse fit than the experimental comparisons (Supplementary Figure 11). Of note, this confirmation of radial dependence in packing argues against other possible alleviators of crypt density such as changes in the size or shape of crypts within mutant clones.

Supplementary Figure 11

Comparison of diffusion model to null model. Results for model fit comparison where the log likelihood of the data under a null model with constant per-clone stromal fraction is compared with the log likelihood of the data under the diffusion model for 1000 potential trajectories (from wild-type to 10-crypt mutant clone) per clone. The black line lies to the left of the density mass, meaning there is evidence for a diffusion-like radial dependence in the crypt packing data.

The inferred diffusion process can be used to define the number of crypt domains impacted to accommodate a new clonal expansion. For example, the model suggests that patches of 10 KDM6A− crypts would require 264 crypt domains to undergo a 1% reduction in their spacing, whereas a 5% reduction would only require 53 crypt domains (Figure 4G and H).

Defining a Homeostatic Threshold

Limited evidence suggests that there are no significant age-related changes in colon length and crypt density. With respect just to STAG2 and KDM6A mutations, the relatively few new crypts arising during life could be easily accommodated by crypt movement. By the time individuals exceed 75 years of age, for every 105 crypts, fission has added only approximately 200 and 290 new STAG2− and KDM6A− crypts, respectively (Figure 4I). However, it seems highly probable that additional genetic variants will also promote fission to different degrees. The potential for diffusion to locally balance this process as fission rate increases was investigated.

Simulations were performed escalating the homeostatic fission rate of 0.7% per year (Supplementary Figure 12). When fission rates remain below approximately 12-fold that of homeostasis, diffusion can generate enough space to accommodate newly generated crypts (Figure 4J). Higher fission rates result in a proportion of clones reaching a threshold of maximum packing density within which crypts are directly touching. This suggests a potential boundary for polyp growth that is dependent on the physical processes of fission and diffusion. We have previously used targeted sequencing of FFPE sections to infer the effect of KRAS activating mutations on crypt fission. Here, in expanding on that initial dataset (see Methods) KRAS activating mutations were found in 35 (22 new) of 256 individuals (130 new) in the age range 20 to 91 years, corresponding to 13.7% of the cohort (Supplementary Figure 13). Mutant allele frequencies in the range of 0.12% to 2.35% combined with total crypt numbers per section enabled estimation of clone size (Supplementary Figure 13). Subsequent mathematical inference indicates that a 17-fold increase in crypt fission rate to 12% (95% CI 10.8–13.7) per year best fits with the data. Approximately 1% of KRAS activating mutations are predicted to breach the threshold for lesion growth after 50 years (Figure 4J). A therapeutic intervention inhibiting crypt fission for any 10 years could reduce this to approximately 0.01% (Figure 4K).

Supplementary Figure 12

Simulating the “breaking point.” Simulations were performed to find the fission rate at which clones may generate new crypts more quickly than can be accommodated by crypt diffusion. Line graphs show the stromal fraction resulting from crypt diffusion and different crypt fission rates in multiples of the homeostatic (wild-type) rate (0.7% per year). Dotted line = stromal fraction calculated from optimal hexagonal packing of circles. r = distance from centroid of patch in crypt domains. Gray area = 95% CI.

Supplementary Figure 13

Evidence for expansion of KRAS-mutant clones in human colon. Targeted amplicon sequencing focusing on KRAS codons 12 and 13 was performed on DNA from FFPE tissue sections from 256 patients of the age range 20 to 91 years. In total, 35 individuals displayed detectable mutations. (A) Mutant allele frequency data of KRAS codons 12 and 13 mutations plotted against age. Data previously published separated from new data obtained as part of this study. (B) KRAS-mutant patch sizes inferred from the mutant allele frequencies and known crypt numbers for 256 individuals plotted against age. (C) Frequency of KRAS mutant crypts for 256 individuals plotted against age. The mean accumulation of mutant clones using the model is plotted in black as well as the 95% CI in gray. Red dotted line shows detection threshold.

Discussion

Several studies have identified large mutant expansions in seemingly normal epithelia.,, In the adult colon, this occurs by increased crypt fission rate, whereby biased mutations can generate large clones that are appropriately distributed within the tissue. In contrast, elevated glandular fission rates are also known to drive the overgrowth of adenomas and CRCs, suggesting that differences in the rate of fission or the response to it must differ between normal and neoplastic tissues.4, 5, 6, 7

In considering the epithelial responses that compensate for elevated fission rates, we first validated a new advantaged clonal mark, KDM6A, that together with STAG2 provided gene-specific assays with 5- and 3-fold increased fission rates, respectively. Comparing the configuration of the size and frequency of clones for both genes demonstrates the different strategies by which age-related mutational burden can be achieved. STAG2 has the higher mutation rate and generates many relatively small clones, whereas KDM6A generates fewer but larger expansions. A corollary of the exponential growth of patches as their size increases is that larger patches will tend to be the most recent and therefore most likely to contain evidence of local adaptation to accommodate new crypts.

The recently recognized process of crypt fusion offers a potential mechanism to compensate for fission. Occurring at equal rates in homeostasis they could effectively balance crypt numbers on a population basis. In considering fusion as a mechanism to accommodate new crypts, a baseline estimate was first established here and found to approximate that for fission; however, no upregulation of fusion accompanying the local expansions resulting from STAG2 and KDM6A mutation was identified. Conceivably other mutations may impact fusion to ease local packing but it does not appear necessary to do so.

Multicrypt clones that form over decades populate the territory previously occupied by multiple independent crypts. Aiming to understand this dispersal, we sought and found evidence of a diffusion-type process in a subset of clones. These are consistent with a recent expansion “caught in the act” of being restored to an ambient crypt density. The behavior captured probably reflects a passive dispersal mechanism rather than actual diffusion and must be accompanied by some level of stromal turnover.

The diffusion coefficient defines the rate of movement of crypt domains and the size of the larger impacted zone that is required to absorb new crypts. Parameterizing the process allows testing of the robustness of the tissue to deal with localized accelerated growth conferred by biased mutations. From this analysis, the homeostatic dispersal mechanism seems able to accommodate increased fission rates of more than 10-fold above baseline. Even for mutations that generate higher crypt fission rates, only the fastest growing clones would overgrow the available space. For example, approximately 5% of clones carrying a gene mutation that confers a 19-fold increase in fission rate would reach a threshold where they lack a stromal domain between crypts and overgrow the available space after 50 years.

The actual threshold at which clonal expansions become recognizable as pathologies may be lower than the extreme one applied here. However, the implication remains that the distinction between phenotypically normal clones and those forming overt pathologies may be determined solely by a probabilistic process in which a recent succession of fission events overwhelms homeostatic dispersal mechanisms.

Activating mutations of KRAS have been described in normal colonic epithelium. The revised estimate of a 17-fold increase for the fission rate conferred by KRAS activating mutations is higher than that previously inferred (10-fold) and is based on analysis of many more patients. It is intriguing that activating mutations of KRAS breach the extreme threshold defined here. KRAS is commonly mutated in a broad spectrum of benign and premalignant pathologies such as serrated lesions and adenomas that may arise at least in part due to the dispersal threshold being reached.29, 30, 31, 32, 33

These findings have clinical significance with respect to bowel cancer screening programs with the implication that clonal expansions with a high malignant potential are not all contained within visible lesions such as sessile serrated adenomas despite having a comparable tissue footprint. Further, it is known that a proportion of adenomas spontaneously regress when observed in longitudinal studies.34, 35, 36 One plausible explanation for such phenomena is that these lesions first develop due to reaching the threshold resulting in local overcrowding but are transient because of ongoing crypt dispersal.

Loss of function mutations affecting the APC tumor suppressor gene are also initiated and expanded by glandular fission.4, 5, 6, 7 Mutation of both APC and KRAS is frequent in CRCs. The 2 pathways are known to interact at the molecular level. It is likely their combined activation will also synergize to further elevate gland fission rate and promote overgrowth as fully neoplastic CRCs develop.

Obesity, a known risk factor for CRC, is known to be accompanied by increased crypt fission rate. Furthermore, diets deficient for methyl donors are known to reduce crypt fission rates in the mouse. An implication of the colon having the capacity to absorb many more new crypts is that modest time limited reductions in fission rates may not only slow the growth of lesions but prevent them from forming at all.