Los sesgos clonales determinan la disponibilidad de mutaciones clave en el cáncer colorrectal para la transformación tumoral.

El estado de las mutaciones impulsoras del cáncer en el tejido normal, específicamente su prevalencia, ventaja selectiva y comportamiento clonal, es fundamental para comprender la etiología del cáncer, ya que las vías disponibles para la transformación dependen de qué mutaciones impulsoras existan y cómo coexistan e interactúen durante la evolución temprana. A pesar de la creciente evidencia de la evolución somática en los tejidos normales, se sabe relativamente poco sobre la carga basal de las mutaciones impulsoras del cáncer en el epitelio colónico humano normal. Lee-Six et al.^[1] informaron que la carga mutacional en el colon está dominada por procesos asociados al envejecimiento y estimaron que aproximadamente el 1% de los criptas en personas de 50 a 60 años contienen al menos una mutación impulsora del cáncer, aunque no se detectaron mutaciones recurrentes en genes canónicos del cáncer colorrectal (CCR), como APC y KRAS. En contraste, estudios dirigidos realizados por nosotros y otros han revelado mutaciones recurrentes en KRAS G12 en el colon normal y han demostrado que la selección positiva influye fuertemente en la carga acumulada de clones mutados en KRAS al promover una expansión clonal sesgada^[2]–^[4]. Sin embargo, aún no está claro si otros genes impulsores del CCR sufren selección en el colon normal, y cuáles son las consecuencias de esto para la adquisición gradual de mutaciones y el orden de los eventos transformadores en la carcinogénesis colorrectal.

Los modelos evolutivos recientes del CCR, incluido el modelo de "Big Bang" del crecimiento tumoral^[5], proponen que muchos de los eventos impulsores clave ocurren temprano, seguidos de una rápida expansión clonal que establece gran parte de la arquitectura intratumoral antes de la detección clínica. Este modelo implica que las etapas más tempranas de la evolución del CCR ya deben involucrar múltiples mutaciones impulsoras, y las estimaciones recientes sugieren que se requieren entre tres y cinco eventos mutacionales para la transformación maligna^[6]–^[9]. Comprender cuáles de estos eventos ya están presentes en el colon normal antes de cualquier cambio histológico o morfológico detectable en el tejido es, por lo tanto, esencial para reconstruir los orígenes del CCR y para distinguir entre la deriva clonal neutral y la evolución impulsada por la selección. La identificación de los eventos tempranos en la evolución somática tiene una relevancia particular para el CCR de inicio temprano, que surge décadas antes de lo esperado^[10] pero que típicamente no muestra evidencia de una carga mutacional marcadamente elevada^[11]. En cambio, el CCR de inicio temprano podría surgir a través de vías evolutivas alternativas en comparación con la enfermedad de inicio tardío^[12]. Además, los fenotipos moleculares de los clones mutados en genes impulsores en el colon normal siguen estando poco caracterizados, y se desconoce si algunos clones ya exhiben características fenotípicas asociadas con un mayor riesgo de transformación. La detección de las firmas moleculares de los clones preneoplásicos puede ofrecer una vía para estratificar las expansiones clonales tempranas según su probabilidad de progresión.

En este trabajo, combinamos el secuenciación dirigida de amplicones, el modelado matemático y la transcriptómica espacial para cuantificar la prevalencia, el comportamiento y el destino evolutivo de los clones mutados en genes impulsores del cáncer en el colon humano normal. Ampliamos el catálogo de mutaciones detectables para incluir variantes en APC, KRAS, TP53, FBXW7 y CTNNB1. Las mutaciones monoalélicas en FBXW7 (R465C) y KRAS (G12) están bajo selección positiva a través de sesgos distintos en la dinámica clonal, mientras que la inactivación monoalélica de APC y TP53 parece ser selectivamente neutral. El modelado predice que la selección altera las trayectorias evolutivas del cáncer, favoreciendo las vías APC-primero en pacientes de mayor edad, al tiempo que permite múltiples vías, incluido KRAS-primero, en individuos más jóvenes. La transcriptómica espacial identifica un subconjunto de clones mutados en KRAS que coexpresan marcadores gástricos e intestinales, lo que indica estados de diferenciación alterados que pueden aumentar el potencial de transformación. En conjunto, estos datos demuestran múltiples vías evolutivas hacia el CCR que surgen en el tejido normal y subrayan la necesidad de estratificar las mutaciones impulsoras tempranas según el riesgo de cáncer.

Resultados

Estrategia para detectar clones mutados en el colon normal

En el epitelio colónico, las mutaciones ventajosas pueden ser seleccionadas positivamente debido a dos procesos independientes^[4],^[13]. En el primero, las células que adquieren mutaciones de novo, junto con su progenie (clones), compiten con las células de criptas de tipo salvaje vecinas para fijarse permanentemente al poblar completamente las glándulas individuales en un proceso conocido como conversión monoclonal^[14]–^[16]. En el segundo, los clones mutados que tienen una ventaja en el proceso de replicación o fisión de las glándulas generan grandes parches multicriptas con el tiempo^[17]–^[20]. Para evaluar las condiciones óptimas para detectar los eventos impulsores del cáncer que promueven cualquiera de los dos procesos, optamos por un análisis de secuenciación basado en amplicones altamente dirigido y replicado de las regiones más frecuentemente mutadas de los genes que se encuentran comúnmente en el CCR (regiones "hotspot" en APC, KRAS, TP53, FBXW7 y CTNNB1) que se pueden aplicar a tejido FFPE en tres escalas diferentes: Micro-seq, AceStrip-seq y Section-seq, lo que garantiza tanto una alta sensibilidad de detección como una amplia cobertura de tejido (Figura suplementaria 1 y "Métodos"). Esto contrasta con el secuenciación del genoma completo de aproximadamente 1000 criptas individuales, adoptado por Lee-Six et al., cuya alta resolución genómica puede identificar más mutaciones por cripta, al tiempo que probablemente omite los eventos fundadores raros de los clones mutados expandidos. El uso de secciones seriadas adyacentes a las utilizadas para el secuenciación de ADN permitió la evaluación histológica y el mapeo espacial de los clones in situ para el análisis posterior (Fig. 1a).

Micro-seq es óptimo para la detección de mutaciones

El análisis de los datos de Micro-seq de 384 muestras (56 pacientes, de 40 a 92 años) que contenían alrededor de 200 criptas microdisecadas por captura láser (N = 76.800 criptas en total) identificó 289 variantes de un solo nucleótido (SNV). AceStrip-seq (43 pacientes, 754.188 criptas) y Section-seq (189 pacientes, 4.114.532 criptas) identificaron 571 y 76 SNV, respectivamente. La comparación del número de SNV detectados por cripta para cada método revela que Micro-seq (3,7 × 10−2 SNV/cripta) es el más eficiente por varios órdenes de magnitud (AceStrip-seq y Section-seq detectaron 1 en 7,5 × 10−4 y 1,9 × 10−5, respectivamente).

Section-seq implica secciones FFPE de área completa que tienen un gran número variable de criptas seccionadas transversalmente, que hemos utilizado previamente para detectar la gran expansión de los clones mutados en KRAS. La baja frecuencia de SNV detectados en todo el panel de amplicones por Section-seq sugiere que tales expansiones grandes son eventos relativamente raros. La menor frecuencia de llamadas de SNV para AceStrip-seq en comparación con Micro-seq indica que, a pesar de utilizar muestras extraídas sistemáticamente de secciones completas para contener alrededor de 1000 criptas por muestra, muchos eventos aún están por debajo del umbral de detección. En consecuencia, los datos de Micro-seq se seleccionaron como la modalidad óptima para la detección de SNV.

Los SNV se acumulan con la edad en el colon normal

Las 384 muestras de Micro-seq, cada una de las cuales contenía alrededor de 200 criptas, se capturaron con láser a partir de secciones FFPE e incluyeron múltiples muestras de cada paciente. El mapeo espacial de las 289 llamadas de mutación de Micro-seq comunes a las submuestras de secciones FFPE se asumió conservadoramente como un solo evento clonal si se derivaban de muestras vecinas, según lo definido por el cálculo de una teselación de Voronoi para cada sección (Fig. 1b, "Métodos" y Tabla suplementaria 1), lo que identificó un total de 285 posibles clones mutados. Las mutaciones en todo el panel de amplicones exhibieron una firma similar a la derivada de las mutaciones codificantes encontradas previamente por el secuenciación del genoma completo de criptas individuales^[1] y dominadas por los cambios CpG > T característicos de la firma similar a un reloj SBS1^[21],^[22] (Fig. 1c, d y Figura suplementaria 2a, b). Aunque se sabe que incluye artefactos asociados con la fijación con formalina, el dominio de SBS1 es una característica esperada de la firma mutacional del epitelio en renovación, incluido el colon^[1]. Además, a pesar del pequeño tamaño de nuestro panel de amplicones (1180 pb), se observó un aumento de los SNV con la edad (5,11 × 10−5 por cripta por año; intervalo de credibilidad del 95% (IC) = 3,8 × 10−5–6,89 × 10−5) (Fig. 1e), lo que proporciona una fuerte evidencia de que la mayoría de las mutaciones son de origen biológico real. Esto equivale a 2,16 mutaciones del exoma por cripta por año (IC del 95% = 1,63–2,91), lo que es comparable, pero ligeramente superior a lo derivado de los datos de Lee-Six et al. (Figura suplementaria 2c, d). Además, de acuerdo con las observaciones previas^[1], Micro-seq reveló una tasa un 20% más alta de acumulación de SNV con la edad en el lado derecho del colon en comparación con el lado izquierdo (Figura suplementaria 2e).

Se encuentran mutaciones de genes impulsores del cáncer en tejido histológicamente normal

Se observaron mutaciones en todos los genes interrogados (Fig. 1f, g). Sin embargo, para CTNNB1, solo se detectaron cuatro mutaciones (incluida G34V, que se observa con frecuencia en el CCR), por lo que este gen se excluyó de los análisis posteriores. Destacaron las mutaciones recurrentes de aminoácidos en las posiciones R465 en FBXW7 y G12 en KRAS (Fig. 1h). Las mutaciones de TP53 se detectaron de forma recurrente en las posiciones R175, R213 y R248, que también son regiones "hotspot" en el CCR. En APC, las mutaciones truncantes recurrentes incluyeron E1353*, Q879* y S1392*, mientras que también se identificaron algunas mutaciones recurrentes (P1409A, R1386T). También se detectó una única mutación truncante R876*, una de las mutaciones APC más comunes en el CCR. Atribuimos la detección de estas mutaciones adicionales en comparación con estudios anteriores^[1] a nuestra estrategia de cribado, que incluye el secuenciación de un número de criptas casi 80 veces mayor.

Las mutaciones detectadas con Micro-seq se compararon con las que se encontraron en los CCR (datos de COSMIC) y se predijo que tenían potencial oncogénico mediante boostDM, un modelo de aprendizaje automático que evalúa los SNV en los genes del cáncer^[23]. Para los cinco genes, se detectó un subconjunto de 45 mutaciones de genes impulsores del cáncer definidas de esta manera en el tejido normal (Fig. 1i). Estos correspondieron a 11 APC, 24 FBXW7, 18 KRAS y 190 clones mutados en genes impulsores del cáncer TP53. Es importante destacar que las mismas y algunas mutaciones adicionales (por ejemplo, APC R1450*) también se detectaron utilizando AceStrip-seq y Section-seq, pero con menor frecuencia (Figura suplementaria 3).

Las regiones de tejido que contenían clones con mutaciones en genes impulsores del cáncer de Micro-seq se revisaron histológicamente. La mayoría se asociaron con tejido histológicamente normal (Fig. 1j). Sin embargo, dos de los 18 clones mutados en KRAS (uno con una mutación G12D y el otro con una mutación G12C) se asociaron con pequeños focos epiteliales con criptas agrandadas, típicas de los pólipos hiperplásicos. Es importante destacar que todos los clones mutados en APC con mutaciones truncantes se encontraron en el epitelio colónico histológicamente normal, lo que destaca que las glándulas mutadas en APC también se pueden encontrar en individuos sin ninguna predisposición genética y que se manifiestan sin ningún cambio morfológico asociado.

Evidencia de selección positiva para las mutaciones de aminoácidos en KRAS y FBXW7

Para evaluar formalmente si los genes mutados están sujetos a selección positiva, se analizó la relación de mutaciones no sinónimas a sinónimas utilizando el método dN/dS adaptado a la región del genoma dirigida en Micro-seq^[24] ("Métodos"). Durante dN/dS, las mutaciones se agrupan para estimar las tasas de mutación de fondo para cada contexto de trinucleótido, lo que también controla los artefactos dentro de una firma de trinucleótido particular. Esto reveló que las mutaciones de aminoácidos en FBXW7 y KRAS están sujetas a una fuerte selección positiva (coeficientes de selección >2 y valores de q ajustados, 3,2 × 10−8 y 1,8 × 10−4, respectivamente). Ningún gen mostró evidencia de selección positiva para las mutaciones truncantes (Fig. 2a). Entre las mutaciones de aminoácidos más frecuentes se encuentran G12D/C/V, Q61L y A146P en KRAS y R465C en FBXW7, lo que sugiere que estos cambios de aminoácidos son los principales impulsores de la selección positiva para las mutaciones en estos genes.

Las mutaciones en FBXW7 y KRAS confieren sesgos a la dinámica clonal

A continuación, analizamos la naturaleza de los sesgos en la dinámica de los clones impulsores del cáncer. Un aumento en la representación tisular de los clones cólicos se logra mediante tasas elevadas de conversión monoclonal (tasa de fijación) o fisión de criptas. Para inferir los valores efectivos de y asociados con grupos de mutaciones, se utilizaron las frecuencias alélicas variantes (VAF) para estimar el número de criptas mutadas contenidas dentro de cada muestra de Micro-seq y se combinaron en un modelo estadístico que considera todos los posibles tamaños de clones mutantes consistentes con los datos observados (Fig. 2b y 'Métodos'). Por ejemplo, en comparación con los grupos de mutaciones sinónimas (20 mutaciones únicas), donde la observación de muchos más eventos contribuye a un valor más alto de , el grupo de mutaciones impulsoras de KRAS G12 (tres mutaciones únicas: G12D, G12C y G12V) presentó un valor inferido de mucho más alto, lo que refleja un sesgo en la fisión de criptas que impulsa la expansión clonal y da como resultado colas más pronunciadas en la distribución predicha del tamaño del clon, lo que corresponde a clones mutantes más grandes en promedio (Fig. 2c, d).

Los valores inferidos de dependen tanto de la tasa de mutación de novo de los contextos trinucleotídicos en los que ocurren, como de las probabilidades conferidas de reemplazo de células madre mutantes (sesgo de fijación). Por lo tanto, las mutaciones dentro de cada grupo se evaluaron individualmente para sus sesgos de fijación, infiriendo valores individuales de y comparándolos entre las mutaciones que surgen en el mismo contexto trinucleotídico, asumiendo que comparten la misma tasa de mutación de novo. FBXW7 R465C fue la única mutación que exhibió un valor individual de sustancialmente más alto que otras mutaciones dentro del mismo contexto trinucleotídico (GCG > T), lo que es consistente con una mayor probabilidad de reemplazo o sesgo de fijación de células madre ( para FBXW7 R465C frente a para una mutación neutra) (Fig. 2e).

La asignación de valores individuales de refinó aún más las estimaciones de la tasa de fisión para cada grupo de mutaciones, confirmando que las tasas elevadas de expansión clonal están asociadas con las mutaciones impulsoras de KRAS G12, como se informó previamente [3,9% de las criptas que sufren fisión por año (IC del 95% = 2,0–7,7%) en comparación con las mutaciones sinónimas (1,1%, IC del 95% = 0,7–1,6%) y las mutaciones impulsoras de KRAS Q61L y A146P combinadas (1,6%, IC del 95% = 0,68–3,1%)]. Se estimó que las mutaciones con sentido de FBXW7 confieren un aumento intermedio en la tasa de fisión (2,5%, IC del 95% = 1,9–3,7%), mientras que los intervalos creíbles de las tasas de fisión para las mutaciones truncantes de APC y las mutaciones impulsoras seleccionadas de TP53 se superpusieron con los de las mutaciones sinónimas (Fig. 2f).

En conjunto, estos resultados sugieren que la selección positiva para las mutaciones impulsoras con sentido de FBXW7 y KRAS se confiere mediante dos modalidades diferentes de dinámica clonal alterada: mayor probabilidad de reemplazo y fijación de células madre para FBXW7 R465C, y tasas muy elevadas de fisión para las mutaciones de KRAS G12.

Los sesgos de selección impactan la carga mutacional

Para estudiar las consecuencias de la selección en la carga mutacional relacionada con la edad en el colon normal, se predijeron las curvas de acumulación basadas en las tasas de fijación y fisión de clones inferidas a partir de los datos de Micro-seq ('Métodos') (Fig. 3). En paralelo, como validación adicional, se aprovecharon y compararon los datos independientes, pero más restringidos, derivados de AceStrip-seq y Section-seq con las predicciones derivadas de Micro-seq (Fig. 3a). Para KRAS y FBXW7, las simulaciones que incluían sesgos en la fisión (KRAS G12) y la fijación (FBXW7 R465C) describieron los datos considerablemente mejor que aquellas sin selección (relación del error cuadrático medio (RRMSE) con/sin selección de 0,68 y 0,42 para KRAS G12 y FBXW7 R465C, respectivamente). No se observó este patrón con las mutaciones impulsoras con sentido de TP53 y, en menor medida, con las mutaciones truncantes de APC (RRMSE 1,03 y 0,82, respectivamente) (Fig. 3a). Por lo tanto, las predicciones realizadas fueron consistentes en las tres escalas de secuenciación.

Luego, comparamos las curvas de carga mutacional entre cada grupo de variantes detectadas, lo que reveló diferencias sorprendentes en la carga mutacional y la prevalencia esperada de diferentes eventos impulsores (Fig. 3b). Por ejemplo, debido a una mayor tasa de fijación, se estima que el 0,1% de todas las criptas contienen la mutación impulsora FBXW7 R465C a la edad media (55 años) en comparación con solo el 0,038% con cualquiera de las mutaciones impulsoras de KRAS G12 (G12D, G12C, G12V), a pesar de que ambos eventos están sujetos a selección positiva. Sin embargo, aunque hay muchos menos eventos mutacionales en general, el alto sesgo de fisión asociado con las mutaciones de KRAS G12 significa que su carga mutacional eventualmente superará la de las mutaciones truncantes de APC o las mutaciones impulsoras con sentido seleccionadas de TP53 a lo largo de la vida (95 años). Por lo tanto, los sesgos de selección dictan la naturaleza y la disponibilidad de las criptas mutantes para la transformación neoplásica como una función no lineal de la edad.

La coocurrencia de mutaciones impulsoras es increíblemente rara en el colon normal

Dado la alta carga mutacional de las mutaciones impulsoras individuales predichas en el colon normal, a continuación, consideramos la coocurrencia de mutaciones, ya que se necesitan múltiples eventos impulsores dentro de la misma cripta para la iniciación del cáncer ^[5]. Cinco de las 384 muestras de Micro-seq compartieron dos instancias de las mutaciones impulsoras estudiadas aquí, pero, en cada caso, las VAF de ambas mutaciones fueron inferiores al 1%, lo que indica que son clones pequeños dentro del área de 200 criptas de LCM recolectada y, por lo tanto, es poco probable que tengan un origen clonal compartido (p < 0,0001, asumiendo que dos o menos criptas albergan alguna de las mutaciones) (Tabla Suplementaria 2). Para explorar la coocurrencia de mutaciones de forma sistemática en áreas más amplias de tejido, optamos por realizar un secuenciación dirigida adicional y la identificación de mutaciones mediante el corte en rodajas de secciones seriadas adyacentes a las muestras de AceStrip-seq donde se detectó una mutación impulsora (Fig. 3c y 'Métodos'). Lógicamente, este enfoque enriquecería los campos de tejido que contienen exclusivamente criptas mutadas, lo que haría que otras mutaciones detectadas en las muestras originales se atribuyeran a eventos de coocurrencia con un alto grado de certeza. Se identificaron y analizaron tres campos grandes de mutaciones impulsoras de KRAS G12V (VAF 22,16, 10,625 y 11,285%), incluido uno de un paciente de 40 años con una carga mutacional excepcionalmente alta. Sin embargo, no se identificó ninguna otra mutación impulsora, lo que sugeriría un posible origen subclonal, dentro de estos campos (Fig. 3c). Por lo tanto, la coocurrencia de dos o más mutaciones impulsoras de la misma cripta parece ser un evento increíblemente raro en el tejido normal, incluso en el contexto de expansiones clonales de KRAS mutantes excepcionalmente grandes.

La selección positiva aumenta el umbral de eventos predicho para la iniciación del CRC

Se ha predicho que el número de eventos mutacionales necesarios para la iniciación del CRC está entre tres y cinco ^[6]–^[9]. Dado que la coocurrencia es rara en el tejido normal, para colocar la dinámica de los clones mutantes observados en este contexto, se adoptó un enfoque basado en simulaciones utilizando un modelo matemático establecido para la iniciación del CRC que involucra las mutaciones impulsoras de APC, KRAS y TP53 ^[25]. El modelo considera todas las posibles rutas a través del espacio de mutación que generan una cripta con inactivación bialélica de los genes supresores de tumores (APC y TP53) y activación monoalélica del oncogén (KRAS) y relaciona la frecuencia de tales criptas con el riesgo de por vida observado de CRC que contienen la misma colección de eventos. Aquí, el modelo se volvió a parametrizar utilizando nuestras tasas de fijación y fisión inferidas para grupos de mutaciones impulsoras monoalélicas detectadas en APC, TP53 y KRAS ('Métodos'). Es importante destacar que existe incertidumbre en el verdadero sesgo de selección de los clones APC nulos, lo que impide una estimación exacta del número de eventos necesarios para la iniciación del CRC. Por lo tanto, para los parámetros restantes no determinados asociados con la inactivación bialélica de APC, consideramos dos posibles escenarios opuestos: ya sea que no se confiera ningún sesgo de selección adicional por un segundo golpe en APC o que tanto los sesgos de fijación como los de fisión estén muy elevados para los clones APC nulos, como se asumió en el modelo "clásico" original (extrapolado a partir de datos de adenomas en ref. ^[25]). Las simulaciones demuestran que, sin un sesgo de selección para los clones APC nulos, la frecuencia de criptas con triple mutación (es decir, aquellas que contienen al menos una mutación detectible en cada uno de los tres genes APC, TP53 y KRAS) se aproxima a la frecuencia de CRC que contienen las mismas mutaciones a la edad de 85 años (Fig. 3d), lo que significa que tres o más eventos serían suficientes para explicar la incidencia de CRC. En contraste, la incorporación de un sesgo de selección para los clones APC nulos da como resultado frecuencias de criptas con triple mutación que exceden la frecuencia de CRC por varios órdenes de magnitud, lo que indica que se requieren eventos adicionales, probablemente alrededor de cinco en total, para la iniciación del CRC (Fig. 3d).

Los eventos fundadores varían con la edad del paciente

A continuación, se utilizaron las simulaciones para explorar el orden temporal de las mutaciones que dan lugar a criptas con triple mutación (criptas que contienen al menos una mutación en cada gen APC, TP53 y KRAS), como posibles precursores del CRC. Basado en los cálculos de carga mutacional anteriores (Fig. 3a) y debido a un gran sesgo de fisión, las variantes de KRAS G12 impulsan las expansiones más grandes de todos los clones que contienen una sola mutación impulsora monoalélica. Sin embargo, las criptas que contienen tres mutaciones impulsoras que coocurren en los tres genes tienen más probabilidades de haber adquirido primero dos golpes independientes en APC a la edad de 80 años, con el 92% de tales clones con triple mutación que surgen con APC como la primera mutación fundadora (Fig. 3e). Como se describió anteriormente ^[25], esta preferencia en el orden temporal se debe a un mayor sesgo de selección asociado con los clones APC nulos que surgen temprano y que confieren una mayor probabilidad de acumular posteriormente múltiples mutaciones impulsoras a lo largo de la vida. Un corolario relevante, dado el aumento de la incidencia de CRC en personas más jóvenes, es que, a una edad más temprana, los clones con triple mutación genética tienen más probabilidades de haber surgido de criptas en las que primero ocurrió una mutación impulsora de KRAS G12 (42,9% frente al 4,7% de los clones a las edades de 50 y 80 años, respectivamente) (Fig. 3e). En conjunto, estas simulaciones sugieren que la selección positiva para la pérdida de APC favorece progresivamente una etiología "APC primero" para el CRC con la edad, pero que, en las personas más jóvenes, son posibles múltiples vías, con KRAS y APC teniendo aproximadamente la misma probabilidad de ser el evento de mutación fundadora en las criptas con triple mutación a la edad de 50 años.

La transcriptómica espacial identifica los clones mutantes de KRAS in situ

La idea de que el CCR podría tener eventos fundacionales alternativos en lugar de proceder necesariamente de una cripta donde la inactivación de APC ocurre primero, nos motivó a analizar con más detalle los clones que portan diferentes mutaciones impulsoras en el tejido normal. Las expansiones clonales impulsadas por mutaciones monoalélicas pueden ser patológicamente indistinguibles del tejido normal, pero poseen diferencias inherentes a nivel molecular que afectan su potencial para acumular mutaciones impulsoras adicionales y su propensión a la progresión neoplásica. Para investigar si los clones con mutaciones impulsoras monoalélicas pueden ser caracterizados fenotípicamente en el tejido histológicamente normal, se realizó una transcriptómica espacial (Visium, 10X) en secciones que eran adyacentes a aquellas que se submuestrearon y que contenían clones mutantes de KRAS [G12V, G12D], TP53 [R175H, R273C], APC [R876*, E1353*] y FBXW7 [R465C] (Figura suplementaria 4). No se observó una agrupación discreta del transcriptoma de las muestras que contenían clones mutantes de APC, TP53 o FBXW7, lo que podría explicarse porque un tamaño de clon relativamente pequeño no genera un cambio transcriptómico suficiente en el conjunto para permitir su detección. Sin embargo, las regiones de tejido que contenían clones mutantes de KRAS G12 formaron grupos distintos del resto de la sección (Fig. 4a, b y Figura suplementaria 5a). Por lo tanto, entre todas las SNV monoalélicas detectadas aquí, solo las mutaciones impulsoras de KRAS G12 se asocian con un fenotipo molecular notablemente distinto del del tejido colónico normal.

La sobreexpresión de REG4 sirve como marcador sustituto para las mutaciones activadoras de KRAS

El distinto grupo fenotípico encontrado en las secciones de tejido que albergan clones mutantes de KRAS G12 mostró una alta expresión del miembro 4 de la familia Regenerating (REG4), lo que permitió la visualización específica de los clones in situ, lo que se confirmó mediante la inmunodetección de la proteína REG4 (Fig. 4c, d y Figura suplementaria 5b). La inmunohistoquímica (IHQ) de REG4 reveló además que el 68% (35/51) de los clones mutantes de KRAS identificados previamente mediante secuenciación de ADN fueron positivos a nivel de proteína, una proporción que aumentó al 92% (22/24) cuando solo se consideraron las mutaciones en la posición de aminoácido G12 (Fig. 4e y Figura suplementaria 5c). Para evaluar la importancia de esta asociación, se realizó una IHQ de REG4 emparejada y una secuenciación de amplicones de KRAS en 16 muestras de tejido que no se habían analizado previamente mediante Micro-seq, identificando diez regiones positivas para REG4 y seis negativas, cada una de las cuales se diseccionó y secuenció mediante microdisección láser para detectar mutaciones de KRAS (Figura suplementaria 5d). Se identificaron mutaciones de KRAS en las posiciones G12 y G13 [G12D, G12C, G12V, G13D] en ocho de las diez (80%) regiones positivas para REG4, mientras que no se identificó ninguna en las seis (0%) regiones negativas para REG4. La prueba exacta de Fisher muestra una asociación altamente significativa (prueba de dos colas, p = 0,007, prueba de una cola, p = 0,0035, que prueba el enriquecimiento de las mutaciones de KRAS en las regiones positivas para REG4) (Fig. 4f). Las regiones positivas para REG4 que no contenían una mutación de KRAS también fueron de tipo salvaje para BRAF V600E, lo cual se sabe que subyace a la neoplasia de pólipos serrados que contienen células con un fenotipo positivo para REG4 comparable ^[26], aunque no se pueden descartar otras mutaciones, como las de NRAS.

Anteriormente, hemos utilizado los tamaños de los clones cuantificados mediante la visualización directa con IHQ para inferir la dinámica de los clones mutantes ^[4]. El uso de la expresión positiva de REG4 como sustituto de la activación de KRAS en el colon normal permite comparar la dinámica de los clones mutantes de KRAS G12 inferida utilizando los tamaños de los clones derivados de VAF con la obtenida mediante la visualización directa utilizando IHQ. Las inferencias basadas en los tamaños de los clones positivos para REG4 indicaron una tasa de fisión elevada del 3,4% por año (IC del 95% = 2,3–4,7%, basado en N = 90 muestras de 53 pacientes) (Fig. 4g), comparable a la inferida para las mutaciones impulsoras de KRAS G12 utilizando los datos de Micro-seq (3,9% por año). En conjunto, la correspondencia en la evaluación cuantitativa de los sesgos de fisión derivados de las mutaciones activadoras de KRAS apoya el uso de la expresión positiva de REG4 como sustituto de la mutación impulsora de KRAS y, de manera más general, la solidez de los datos de secuenciación utilizados aquí para identificar dichos sesgos de selección.

Los clones mutantes de KRAS coexpresan marcadores de linaje intestinal y gástrico

Las predicciones del modelo de que las mutaciones impulsoras de KRAS G12 son eventos fundacionales candidatos para las criptas con triple mutación, especialmente en individuos más jóvenes, nos motivaron a investigar las características moleculares de los clones mutantes de KRAS e identificar si se podían detectar elementos que recordaran al CCR en esta etapa temprana. Utilizando la expresión positiva de REG4 para identificar 14 clones mutantes de KRAS capturados por el perfilado de Visium, identificamos una firma mutacional de KRAS que comprende 628 genes. Entre los genes expresados diferencialmente se encontraban MUC2, TFF1, S100P y SPINK4 (Fig. 4h). El análisis de enriquecimiento de conjuntos de genes (GSEA) mostró un enriquecimiento de tipos de células epiteliales no intestinales, incluidos los tipos de células esofágicas, gástricas y fetales (Fig. 4i). Las firmas transcriptómicas asociadas con estados celulares que expresan marcadores fetales o metaplásicos ^[27]–^[30] también se enriquecieron dentro de la firma de clones mutantes de KRAS (Figura suplementaria 5e). En particular, no hubo enriquecimiento de las firmas clásicas de células madre intestinales ^[30] y linajes secretores del colon ^[31], ^[32].

La metaplasia gástrica se ha informado previamente en pólipos serrados con mutaciones de BRAF definidas por la sobreexpresión de MUC5AC y la pérdida del marcador intestinal CDX2 ^[26]. Los clones mutantes de KRAS expresaron marcadores gástricos, incluido MUC5AC (Figura suplementaria 5f), pero la conversión completa del linaje gástrico (MUC5AC+/CDX2−) solo se observó en un clon, y el 60% de los clones mutantes de KRAS coexpresaron MUC5AC y CDX2 (Fig. 4j). El tinte de inmunofluorescencia de ocho clones mutantes de KRAS que coexpresan MUC5AC/CDX2 destacó que, en todos los clones, alrededor del 5% de las células mutantes de KRAS expresaban simultáneamente marcadores gástricos e intestinales (Fig. 4k). En el estómago y el esófago, se ha descrito una coexpresión comparable de marcadores específicos de tejidos como metaplasia incompleta o confusión de linaje, y se ha asociado con un alto riesgo de progresión a displasia y cáncer ^[33]. Aunque la mezcla de linajes observada aquí podría ser un estado transitorio, hasta donde sabemos, este fenómeno no se ha descrito previamente en el epitelio colónico normal.

Finalmente, analizamos la relevancia de esta firma de KRAS mutante para el CCR. La metaplasia gástrica se ha relacionado con el subtipo transcripcional iCMS3 del CCR, donde REG4 y SPINK4 son genes de firma definitorios de esta clasificación, y entre los cuales existe un enriquecimiento conocido de tumores mutados en KRAS ^[34]. Aquí, observamos un enriquecimiento de los marcadores CMS3/iCMS3 en los clones mutantes de KRAS (Fig. 4l), lo que sugiere que estos perfiles iCMS3 pueden ser una consecuencia de un evento de mutación impulsora temprana de KRAS que confiere características gástricas que se mantienen a lo largo de la evolución de este subconjunto de CCR.

Estrategia para detectar clones mutantes en el colon normal

En el epitelio colónico, las mutaciones ventajosas pueden ser seleccionadas positivamente debido a dos procesos independientes ^[4], ^[13]. En el primero, las células que adquieren mutaciones de novo, junto con su progenie posterior (clones), compiten con las células de tipo salvaje vecinas de la cripta para convertirse en fijas permanentemente al poblar completamente las glándulas individuales en un proceso conocido como conversión monoclonal ^[14]–^[16]. En el segundo, los clones mutantes con ventaja en el proceso de replicación o fisión de las glándulas generan grandes parches multicriptas con el tiempo ^[17]–^[20]. Para evaluar las condiciones óptimas para detectar eventos de mutación que promuevan cualquiera de los dos procesos, optamos por un análisis de secuenciación basado en réplicas altamente dirigido de las regiones más frecuentemente mutadas de los genes que se encuentran comúnmente en el CCR (regiones de puntos calientes en APC, KRAS, TP53, FBXW7 y CTNNB1) que se podrían aplicar a tejido FFPE en tres escalas diferentes, Micro-seq, AceStrip-seq y Section-seq, lo que garantiza tanto una alta sensibilidad de detección como una amplia cobertura de tejido (Figura suplementaria 1 y "Métodos"). Esto contrasta con la secuenciación del genoma completo de aproximadamente 1000 criptas individuales adoptada por Lee-Six et al., cuya alta resolución genómica puede identificar más mutaciones por cripta, pero probablemente omite eventos raros y fundacionales de clones mutantes expandidos. El uso de secciones seriadas adyacentes a las utilizadas para la secuenciación de ADN permitió la evaluación histológica y el mapeo espacial de los clones in situ para el análisis posterior (Fig. 1a).

Micro-seq es óptimo para la detección de mutaciones

El análisis de los datos de Micro-seq de 384 muestras (56 pacientes, 40–92 años), cada una de las cuales contenía alrededor de 200 criptas microdisecadas mediante captura láser (N = 76 800 criptas en total), identificó 289 variantes de nucleótido único (SNV). AceStrip-seq (43 pacientes, 754 188 criptas) y Section-seq (189 pacientes, 4 114 532 criptas) identificaron 571 y 76 SNV, respectivamente. La comparación del número de SNV detectados por cripta para cada método revela que Micro-seq (3,7 × 10−2 SNV/cripta) es el más eficiente por varios órdenes de magnitud (AceStrip-seq y Section-seq detectaron 1 en 7,5 × 10−4 y 1,9 × 10−5, respectivamente).

Section-seq implica secciones FFPE de área completa que tienen un gran número variable de criptas seccionadas, que hemos utilizado previamente para detectar una gran expansión de clones mutantes de KRAS. La baja frecuencia de SNV detectados en todo el panel de amplicones por Section-seq sugiere que tales expansiones grandes son eventos relativamente raros. La menor frecuencia de llamadas de SNV para AceStrip-seq en comparación con Micro-seq indica que, a pesar de utilizar muestras extraídas sistemáticamente de secciones completas para contener alrededor de 1000 criptas por muestra, muchos eventos aún están por debajo del umbral de detección. En consecuencia, los datos de Micro-seq se seleccionaron como la modalidad óptima para la detección de SNV.

Las SNV se acumulan con la edad en el colon normal

Las 384 muestras de Micro-seq, cada una de las cuales contenía alrededor de 200 criptas, se capturaron mediante láser a partir de secciones FFPE e incluyeron múltiples muestras de cada paciente. El mapeo espacial de las 289 llamadas de mutación de Micro-seq comunes a las submuestras de secciones FFPE se asumió de forma conservadora como un único evento clonal si se derivaban de muestras vecinas, según lo definido por el cálculo de una teselación de Voronoi para cada sección (Fig. 1b, "Métodos" y Tabla suplementaria 1), lo que identificó un total de 285 clones mutantes putativos. Las mutaciones en todo el panel de amplicones exhibieron una firma similar a la derivada de las mutaciones codificantes encontradas previamente mediante la secuenciación del genoma completo de criptas individuales ^[1] y dominada por los cambios de CpG > T característicos de la firma similar a un reloj SBS1 ^[21], ^[22] (Fig. 1c, d y Figura suplementaria 2a, b). Aunque se sabe que incluye artefactos asociados con la fijación con formalina, el dominio de SBS1 es una característica esperada de la firma de mutación del epitelio en renovación, incluido el colon ^[1]. Además, a pesar del pequeño tamaño de nuestro panel de amplicones (1180 pb), se observó un aumento de las SNV con la edad (5,11 × 10−5 por cripta por año; IC del 95% = 3,8 × 10−5–6,89 × 10−5) (Fig. 1e), lo que proporciona una fuerte evidencia de que la mayoría de las mutaciones son de origen biológico real. Esto equivale a 2,16 mutaciones de exoma por cripta por año (IC del 95% = 1,63–2,91), lo cual es comparable, pero ligeramente superior al derivado de los datos de Lee-Six et al. (Figura suplementaria 2c, d). Además, de acuerdo con las observaciones anteriores ^[1], Micro-seq reveló una tasa un 20% mayor de acumulación de SNV con la edad en el lado derecho del colon en comparación con el izquierdo (Figura suplementaria 2e).

Las mutaciones impulsoras del cáncer se encuentran en el tejido histológicamente normal

Se observaron mutaciones en todos los genes analizados (Fig. 1f, g). Sin embargo, para CTNNB1, solo se detectaron cuatro mutaciones (incluida G34V, que se observa con frecuencia en el CCRC), por lo que este gen se excluyó de los análisis posteriores. Destacaron las mutaciones recurrentes con cambio de aminoácido en las posiciones R465 de FBXW7 y G12 de KRAS (Fig. 1h). Las mutaciones de TP53 se detectaron de forma recurrente en las posiciones R175, R213 y R248, que también son regiones de "hotspot" en el CCRC. En APC, las mutaciones truncantes recurrentes incluyeron E1353*, Q879* y S1392*, mientras que también se identificaron algunas mutaciones recurrentes con cambio de aminoácido (P1409A, R1386T). También se detectó una única mutación truncante R876*, una de las mutaciones más comunes de APC en el CCRC. Atribuimos la detección de estas mutaciones adicionales en comparación con estudios previos^[1] a nuestra estrategia de cribado, que incluye el secuenciamiento de un número de criptas aproximadamente 80 veces mayor.

Las mutaciones detectadas con Micro-seq se compararon con las encontradas en los CCRC (datos de COSMIC) y se predijo que tendrían potencial oncogénico mediante boostDM, un modelo de aprendizaje automático que evalúa las SNV en los genes del cáncer^[23]. Para los cinco genes, se detectó un subconjunto de 45 mutaciones de genes impulsores del cáncer de CRC en el tejido normal (Fig. 1i). Estas correspondieron a 11 APC, 24 FBXW7, 18 KRAS y 190 clones mutantes del gen TP53 que actúan como impulsores del cáncer. Es importante destacar que las mismas y algunas mutaciones adicionales (por ejemplo, APC R1450*) que actúan como impulsores del cáncer también se detectaron utilizando AceStrip-seq y Section-seq, pero con frecuencias más bajas (Fig. S3).

Se revisaron histológicamente las regiones de tejido que contenían clones con mutaciones de genes impulsores del cáncer detectadas con Micro-seq. La mayoría se asociaron con tejido histológicamente normal (Fig. 1j). Sin embargo, dos de los 18 clones mutantes de KRAS (uno con una mutación G12D y el otro con una mutación G12C) se asociaron con pequeños focos epiteliales con criptas agrandadas, típicas de los pólipos hiperplásicos. Es importante destacar que todos los clones mutantes de APC con mutaciones truncantes se encontraron en el epitelio colónico histológicamente normal, lo que pone de manifiesto que los glándulas mutantes de APC también se pueden encontrar en individuos sin ninguna predisposición genética y manifestarse sin cambios morfológicos asociados.

Evidencia de selección positiva para las mutaciones con cambio de aminoácido de KRAS y FBXW7

Para evaluar formalmente si los genes mutados están sujetos a selección positiva, se analizó la relación entre las mutaciones no sinónimas y las sinónimas utilizando el método dN/dS adaptado a la región del genoma analizada en Micro-seq^[24] ("Métodos"). Durante el análisis dN/dS, las mutaciones se agrupan para estimar las tasas de mutación de fondo para cada contexto de trinucleótido, lo que también controla los artefactos dentro de una firma de trinucleótido particular. Esto reveló que las mutaciones con cambio de aminoácido en FBXW7 y KRAS están sujetas a una fuerte selección positiva (coeficientes de selección >2 y valores q ajustados, 3,2 × 10−8 y 1,8 × 10−4, respectivamente). Ningún gen mostró evidencia de selección positiva para las mutaciones truncantes (Fig. 2a). Entre las mutaciones con cambio de aminoácido más frecuentes se encuentran G12D/C/V, Q61L y A146P en KRAS y R465C en FBXW7, lo que sugiere que estos cambios de aminoácido son los principales impulsores de la selección positiva para las mutaciones en estos genes.

Las mutaciones de FBXW7 y KRAS confieren sesgos a la dinámica de los clones

A continuación, buscamos la naturaleza de los sesgos en la dinámica de los clones impulsores del cáncer. Un aumento en la representación del tejido de los clones colónicos se logra mediante tasas elevadas de conversión monoclonal (tasa de fijación) o fisión de criptas (). Para inferir los valores efectivos de y asociados con grupos de mutaciones, se utilizaron las frecuencias alélicas variantes (VAF) para estimar el número de criptas mutadas contenidas dentro de cada muestra de Micro-seq y se combinaron en un modelo estadístico que considera todos los posibles tamaños de clones mutantes consistentes con los datos observados (Fig. 2b y "Métodos"). Por ejemplo, en comparación con las mutaciones sinónimas agrupadas (20 mutaciones únicas), donde la observación de muchos más eventos contribuye a un valor más alto de , el grupo de mutaciones del gen KRAS G12 (tres mutaciones únicas: G12D, G12C y G12V) presentó un valor inferido de mucho más alto, lo que refleja un sesgo en la fisión de criptas que impulsa la expansión clonal y da como resultado colas más pesadas en la distribución del tamaño de los clones predichos, lo que corresponde a clones mutantes más grandes en promedio (Fig. 2c, d).

Los valores inferidos de dependen tanto de la tasa de mutación de novo de los contextos de trinucleótido en los que ocurren como de las probabilidades conferidas de reemplazo de las células madre mutantes (sesgo de fijación). Por lo tanto, las mutaciones dentro de cada grupo se evaluaron individualmente para sus sesgos de fijación infiriendo valores individuales de y comparándolos entre las mutaciones que surgen en el mismo contexto de trinucleótido, asumiendo que comparten la misma tasa de mutación de novo. FBXW7 R465C fue la única mutación que exhibió un valor individual de sustancialmente más alto que otras mutaciones dentro del mismo contexto de trinucleótido (GCG > T), lo que es consistente con una mayor probabilidad de reemplazo o sesgo de fijación de las células madre ( para FBXW7 R465C frente a para una mutación neutra) (Fig. 2e).

La asignación de valores individuales de refinó aún más las estimaciones de la tasa de fisión para cada grupo de mutaciones, lo que confirmó que las tasas elevadas de expansión clonal están asociadas con las mutaciones del gen KRAS G12, como se informó previamente [3,9% de las criptas que sufren fisión por año (IC del 95% = 2,0-7,7%) en comparación con las mutaciones sinónimas (1,1%, IC del 95% = 0,7-1,6%) y las mutaciones del gen KRAS Q61L y A146P combinadas (1,6%, IC del 95% = 0,68-3,1%)]. Las mutaciones con cambio de aminoácido de FBXW7 se estimaron para conferir un aumento intermedio en la tasa de fisión (2,5%, IC del 95% = 1,9-3,7%), mientras que los intervalos de credibilidad de las tasas de fisión para las mutaciones truncantes de APC y las mutaciones seleccionadas del gen TP53 se superpusieron con los de las mutaciones sinónimas (Fig. 2f).

En conjunto, estos resultados sugieren que la selección positiva para las mutaciones con cambio de aminoácido de FBXW7 y KRAS se confiere mediante dos modalidades diferentes de alteración de la dinámica de los clones: aumento de la probabilidad de reemplazo de las células madre y fijación para FBXW7 R465C, y tasas muy elevadas de fisión para las mutaciones del gen KRAS G12.

Los sesgos de selección impactan la carga mutacional

Para estudiar las consecuencias de la selección en el colon normal, se predijeron las curvas de acumulación de la carga mutacional relacionada con la edad basándose en las tasas de fijación y fisión de los clones inferidas a partir de los datos de Micro-seq ("Métodos") (Fig. 3). En paralelo, como validación adicional, se aprovecharon los datos independientes pero más restringidos derivados de AceStrip-seq y Section-seq y se compararon con las predicciones derivadas de Micro-seq (Fig. 3a). Para KRAS y FBXW7, las simulaciones que incluyen sesgos en la fisión (KRAS G12) y la fijación (FBXW7 R465C) funcionaron considerablemente mejor para describir los datos que aquellas sin selección (relación del error cuadrático medio (RRMSE) con/sin selección 0,68 y 0,42 para KRAS G12 y FBXW7 R465C, respectivamente). Este patrón no se observó con las mutaciones con cambio de aminoácido del gen TP53 y, en menor medida, con las mutaciones truncantes de APC (RRMSE 1,03 y 0,82, respectivamente) (Fig. 3a). Por lo tanto, las predicciones realizadas fueron consistentes en las tres escalas de secuenciación.

A continuación, se compararon las curvas de carga mutacional entre cada grupo de variantes detectadas, lo que reveló diferencias sorprendentes en la carga mutacional y la prevalencia esperada de diferentes eventos impulsores (Fig. 3b). Por ejemplo, debido a una mayor tasa de fijación, se estima que el 0,1% de todas las criptas contienen la mutación del gen FBXW7 R465C a la edad media (55 años) en comparación con solo el 0,038% con cualquiera de las mutaciones del gen KRAS G12 (G12D, G12C, G12V), a pesar de que ambos eventos están sujetos a selección positiva. Sin embargo, aunque hay muchos menos eventos mutacionales en general, el alto sesgo de fisión asociado con las mutaciones del gen KRAS G12 significa que su carga mutacional eventualmente superará la de las mutaciones truncantes detectadas de APC o las mutaciones seleccionadas del gen TP53 a lo largo de la vida (95 años) (Fig. 3b). Los sesgos de selección dictan, por lo tanto, la naturaleza y la disponibilidad de los clones mutantes para la transformación neoplásica como una función no lineal de la edad.

La coocurrencia de mutaciones impulsoras es increíblemente rara en el colon normal

Dado el alto número de mutaciones de genes impulsores individuales predichas en el colon normal, a continuación, se consideró la coocurrencia de mutaciones, ya que se necesitan múltiples eventos impulsores dentro del mismo cripto para la iniciación del cáncer^[5]. Cinco de las 384 muestras de Micro-seq compartieron dos instancias de las mutaciones del gen impulsor estudiadas aquí, pero, en cada caso, las VAF de ambas mutaciones fueron inferiores al 1%, lo que indica que son clones pequeños dentro del área de 200 criptas obtenida por disección láser por captura de microdiseción (LCM), y por lo tanto es poco probable que tengan un origen clonal compartido (p < 0,0001, asumiendo que dos o menos criptas albergan cualquiera de las mutaciones) (Tabla S2). Para explorar sistemáticamente la coocurrencia de genes impulsores en porciones más amplias de tejido, se optó por realizar un secuenciamiento dirigido adicional y la identificación de mutaciones mediante el corte en rodajas de secciones seriadas adyacentes a las muestras de AceStrip-seq donde se detectó una mutación del gen impulsor (Fig. 3c y "Métodos"). Lógicamente, este enfoque enriquecería los campos de tejido que contienen exclusivamente criptas mutadas, lo que haría que otras mutaciones del gen impulsor detectadas en las muestras originales fueran atribuibles a eventos de coocurrencia con un alto grado de certeza. Se identificaron y analizaron tres grandes campos de mutaciones del gen KRAS G12V (VAF 22,16, 10,625 y 11,285%), incluido uno de un paciente de 40 años con una carga mutacional excepcionalmente alta. Sin embargo, no se identificó ninguna otra mutación del gen impulsor, lo que sugeriría un posible origen subclonal (Fig. 3c). Por lo tanto, la coocurrencia de dos o más mutaciones del gen impulsor del mismo cripto parece ser un evento increíblemente raro en el tejido normal, incluso en el contexto de expansiones clonales de KRAS muy grandes.

La selección positiva aumenta el umbral de eventos predicho para la iniciación del CCRC

Se ha predicho que el número de eventos mutacionales necesarios para la iniciación del CCRC se sitúa entre tres y cinco^[6]-^[9]. Dado que la coocurrencia es rara en el tejido normal, para contextualizar la dinámica de los clones mutantes observados, se adoptó un enfoque basado en simulaciones utilizando un modelo matemático establecido para la iniciación del CCRC que involucra las mutaciones de los genes APC, KRAS y TP53^[25]. El modelo considera todas las posibles rutas a través del espacio de mutación que generan un cripto con inactivación bialélica de los genes supresores de tumores (APC y TP53) y activación monoalélica del oncogén (KRAS) y relaciona la frecuencia de dichos criptos con el riesgo de por vida observado de los CCRC que contienen la misma colección de eventos. Aquí, el modelo se volvió a parametrizar utilizando nuestras tasas de fijación y fisión inferidas para grupos de mutaciones monoalélicas detectadas en APC, TP53 y KRAS ("Métodos"). Es importante destacar que existe incertidumbre en el verdadero sesgo de selección de los clones APC nulos, lo que impide una estimación exacta del número de eventos necesarios para la iniciación del CCRC. Por lo tanto, para los parámetros restantes no determinados asociados con la inactivación bialélica de APC, se consideraron dos posibles escenarios opuestos: ya sea que no se confiera ningún sesgo de selección adicional por un segundo evento en APC o que tanto los sesgos de fijación como los de fisión se eleven en gran medida para los clones APC nulos, como se asumió en el modelo "clásico" original (extrapolado a partir de los datos de adenomas en la ref. ^[25]). Las simulaciones demuestran que, sin un sesgo de selección para los clones APC nulos, la frecuencia de los criptos con triple mutación (es decir, aquellos que contienen al menos una mutación del gen impulsor detectable en cada uno de los tres genes APC, TP53 y KRAS) se aproxima a la frecuencia de los CCRC que contienen las mismas mutaciones a la edad de 85 años (Fig. 3d), lo que significa que tres o más eventos serían suficientes para explicar la incidencia del CCRC. En contraste, la incorporación de un sesgo de selección para los clones APC nulos da como resultado frecuencias de criptos con triple mutación que exceden la frecuencia de los CCRC en varios órdenes de magnitud, lo que indica que se requieren eventos adicionales, probablemente alrededor de cinco en total, para la iniciación del CCRC (Fig. 3d).

Los eventos fundadores varían con la edad del paciente

A continuación, se utilizaron las simulaciones para explorar el orden temporal de las mutaciones que dan lugar a criptas con triple mutación (criptas que contienen al menos una mutación en cada gen APC, TP53 y KRAS), como posibles precursoras del CRC. Basándose en cálculos previos de la carga mutacional (Fig. 3a) y debido a un gran sesgo de fisión, las variantes de KRAS G12 impulsan las mayores expansiones de todos los clones que contienen una única mutación monoalélica. Sin embargo, las criptas que contienen tres mutaciones conductoras que ocurren simultáneamente en los tres genes tienen más probabilidades de haber adquirido primero dos mutaciones independientes en APC a la edad de 80 años, y el 92% de estos clones con triple mutación surgen con APC como la primera mutación fundadora (Fig. 3e). Como se describió anteriormente^[25], esta preferencia en el orden temporal se debe a un mayor sesgo de selección asociado con los clones APC nulos que surgen al principio y que confieren una mayor probabilidad de acumular posteriormente múltiples mutaciones conductoras a lo largo de la vida. Una implicación relevante, dado el aumento de la incidencia de los CRC en individuos más jóvenes, es que, a una edad más temprana, los clones con mutaciones en tres genes tienen más probabilidades de haber surgido de criptas en las que primero se produjo una mutación conductora KRAS G12 (42,9% frente al 4,7% de los clones a las edades de 50 y 80 años, respectivamente) (Fig. 3e). En conjunto, estas simulaciones sugieren que la selección positiva para la pérdida de APC favorece progresivamente una etiología "APC primero" para el CRC con la edad, pero que, en los individuos más jóvenes, son posibles múltiples vías, siendo KRAS y APC aproximadamente igualmente probables como el evento de mutación conductora fundadora en las criptas con triple mutación a la edad de 50 años.

La transcriptómica espacial identifica clones con mutación KRAS in situ

La idea de que el CRC podría tener eventos fundadores alternativos en lugar de proceder necesariamente de una cripta en la que primero se produjo la inactivación de APC, nos motivó a examinar con más detalle los clones que portan diferentes mutaciones conductoras en el tejido normal. Las expansiones clonales impulsadas por mutaciones monoalélicas pueden ser patológicamente indistinguibles del tejido normal, pero poseen diferencias inherentes a nivel molecular que afectan a su potencial para acumular mutaciones conductoras adicionales y a su propensión a la progresión neoplásica. Para investigar si los clones con mutaciones conductoras monoalélicas pueden ser perfilados fenotípicamente en el tejido histológicamente normal, se realizó una transcriptómica espacial (Visium, 10X) en secciones que eran adyacentes a las que se submuestrearon y que contenían clones con mutación KRAS [G12V, G12D], TP53 [R175H, R273C], APC [R876*, E1353*] y FBXW7 [R465C] (Fig. S4). No se observó una agrupación discreta del transcriptoma de las muestras que contenían clones con mutación APC, TP53 o FBXW7, lo que podría explicarse porque un tamaño de clon relativamente pequeño no genera un cambio transcriptómico suficiente en relación con el volumen total para permitir su detección. Sin embargo, las regiones de tejido que contenían clones con mutación conductora KRAS G12 formaron grupos distintos del resto de la sección (Fig. 4a, b y Fig. S5a). Por lo tanto, entre todas las SNV monoalélicas detectadas aquí, solo las mutaciones conductoras KRAS G12 están asociadas con un fenotipo molecular notablemente distinto del del tejido colónico normal.

La sobreexpresión de REG4 sirve como marcador sustituto de las mutaciones activadoras de KRAS

El distinto grupo fenotípico que se encontró en las secciones de tejido que albergaban clones con mutación KRAS G12 mostró una alta expresión de Regenerating family member 4 (REG4), lo que permitió la visualización específica de los clones in situ, lo que se confirmó mediante la inmunodetección de la proteína REG4 (Fig. 4c, d y Fig. S5b). La inmunohistoquímica (IHC) de REG4 reveló además que el 68% (35/51) de los clones con mutación KRAS identificados previamente mediante secuenciación de ADN eran positivos a nivel de proteína, proporción que aumentó al 92% (22/24) cuando solo se consideraron las mutaciones en la posición de aminoácido G12 (Fig. 4e y Fig. S5c). Para evaluar la importancia de esta asociación, se realizó una IHC de REG4 emparejada con la secuenciación de amplicones de KRAS en 16 muestras de tejido que no se habían analizado previamente mediante Micro-seq, identificando diez regiones positivas para REG4 y seis negativas, cada una de las cuales se diseccionó y secuenció mediante microdisección láser para detectar mutaciones de KRAS (Fig. S5d). Se identificaron mutaciones de KRAS en las posiciones G12 y G13 [G12D, G12C, G12V, G13D] en ocho de las diez (80%) regiones positivas para REG4, mientras que no se identificó ninguna en las seis (0%) regiones negativas para REG4. La prueba exacta de Fisher muestra una asociación altamente significativa (prueba de dos colas p = 0,007, prueba de una cola p = 0,0035, que prueba el enriquecimiento de las mutaciones de KRAS en las regiones positivas para REG4) (Fig. 4f). Las regiones positivas para REG4 que no contenían una mutación de KRAS también eran de tipo salvaje para BRAF V600E, que también se sabe que subyace a la neoplasia de pólipos serrados que contienen células con un fenotipo comparable positivo para REG4^[26], aunque no se pueden descartar otras mutaciones, como las de NRAS.

Anteriormente, hemos utilizado tamaños de clon cuantificados mediante la visualización directa con IHC para inferir la dinámica de los clones mutados^[4]. El uso de la expresión positiva de REG4 como sustituto de la activación de KRAS en el colon normal permite comparar la dinámica de los clones con mutación KRAS G12 inferida mediante tamaños de clon derivados de VAF con la obtenida mediante la visualización directa con IHC. Las inferencias basadas en los tamaños de los clones positivos para REG4 indicaron una tasa de fisión elevada del 3,4% por año (IC del 95% = 2,3-4,7%, basada en N = 90 muestras de 53 pacientes) (Fig. 4g), comparable a la inferida para las mutaciones conductoras KRAS G12 utilizando datos de Micro-seq (3,9% por año). En conjunto, la correspondencia en la evaluación cuantitativa de los sesgos de fisión derivados de las mutaciones activadoras de KRAS apoya tanto la expresión positiva de REG4 como sustituto de la mutación conductora KRAS y, en general, la solidez de los datos de secuenciación utilizados aquí para identificar dichos sesgos de selección.

Los clones con mutación KRAS coexpresan marcadores de linaje intestinal y gástrico

Las predicciones del modelo de que las mutaciones conductoras KRAS G12 son eventos fundadores candidatos para las criptas con triple mutación, especialmente en los individuos más jóvenes, nos motivaron a investigar las características moleculares de los clones con mutación KRAS e identificar si se podían detectar elementos que recordaran al CRC en esta etapa temprana. Utilizando la expresión positiva de REG4 para identificar 14 clones con mutación KRAS capturados mediante el perfilado de Visium, identificamos una firma de mutación KRAS que comprende 628 genes. Entre los genes expresados diferencialmente se encontraban MUC2, TFF1, S100P y SPINK4 (Fig. 4h). El análisis de enriquecimiento de conjuntos de genes (GSEA) mostró un enriquecimiento de tipos de células epiteliales no intestinales, incluidos los tipos de células esofágicas, gástricas y fetales (Fig. 4i). Las firmas transcriptómicas asociadas con estados celulares que expresan marcadores fetales o metaplásicos^[27]–^[30] también se enriquecieron dentro de la firma de clones con mutación KRAS (Fig. S5e). En particular, no hubo enriquecimiento de las firmas clásicas de células madre intestinales^[30] y de los linajes secretores del colon^[31],^[32].

La metaplasia gástrica se ha informado previamente en los pólipos serrados con mutaciones BRAF definidas por la sobreexpresión de MUC5AC y la pérdida del marcador intestinal CDX2^[26]. Los clones con mutación KRAS expresaron marcadores gástricos, incluido MUC5AC (Fig. S5f), pero la conversión completa del linaje gástrico (MUC5AC+/CDX2−) solo se observó en un clon, y el 60% de los clones con mutación KRAS coexpresaron MUC5AC y CDX2 (Fig. 4j). El tinte de inmunofluorescencia de ocho clones con mutación KRAS que coexpresan MUC5AC/CDX2 destacó que, en todos los clones, alrededor del 5% de las células con mutación KRAS expresaban simultáneamente marcadores gástricos e intestinales (Fig. 4k). En el estómago y el esófago, se ha descrito una coexpresión comparable de marcadores específicos de tejidos como metaplasia incompleta o confusión de linajes, y se ha asociado con un alto riesgo de progresión a displasia y cáncer^[33]. Aunque la mezcla de linajes que se observa aquí podría ser un estado transitorio, hasta donde sabemos, este fenómeno no se ha descrito previamente en el epitelio colónico normal.

Finalmente, analizamos la relevancia de esta firma de KRAS mutado para el CRC. La metaplasia gástrica se ha relacionado con el subtipo transcripcional iCMS3 del CRC, donde REG4 y SPINK4 son genes de firma definitorios de esta clasificación, y entre los cuales existe un enriquecimiento conocido de tumores con mutación KRAS^[34]. Aquí, observamos un enriquecimiento de los marcadores CMS3/iCMS3 en los clones con mutación KRAS (Fig. 4l), lo que sugiere que estos perfiles iCMS3 pueden ser una consecuencia de un evento de mutación conductora KRAS temprana que confiere características gástricas que se mantienen a lo largo de la evolución de este subconjunto de CRC.

Discusión

Utilizando un enfoque de secuenciación de amplicones basado en réplicas altamente dirigido, la lista de mutaciones conductoras del cáncer que se encuentran en el colon normal se ha ampliado para incluir mutaciones clave del CRC, incluidas R876*/R1450* en APC, R465C en FBXW7, G12C/D/V en KRAS y R175H/R248Q/W en TP53. En particular, se encontraron mutaciones truncantes de APC en el tejido histológicamente normal en individuos sin predisposición genética al CRC familiar.

Se descubrió que las mutaciones monoalélicas en KRAS y FBXW7, pero no en APC o TP53, mostraban una fuerte selección positiva en el epitelio colónico normal mediada por diferentes comportamientos clonales. Las mutaciones FBXW7 R645C se ven favorecidas por la competencia de células madre sesgada, lo que conduce a tasas más altas de fijación dentro de la cripta. Los clones que portan mutaciones activadoras de KRAS en la posición de aminoácido G12 se expanden preferentemente a través de tasas más altas de fisión de la cripta. Anteriormente, en ratones, se ha demostrado que las mutaciones KRAS tienen sesgos tanto en la fijación como en la expansión de los clones^[19],^[35],^[36]. Aquí, no encontramos evidencia clara de lo último en el epitelio humano, lo que sugiere una posible diferencia de especies en el impacto, pero los datos comparables en humanos siguen siendo limitados. Dado su sesgo de selección en el tejido normal, tanto las mutaciones activadoras de KRAS como las mutaciones FBXW7 R465C tienen una mayor biodisponibilidad como eventos tempranos en el CRC, mientras que los resultados neutrales asociados con las mutaciones truncantes de APC de un solo impacto y las mutaciones de TP53 sugieren que su carga mutacional acumulativa está dictada principalmente por su susceptibilidad a los procesos mutacionales activos en el epitelio colónico. Por lo tanto, APC requiere una inactivación bialélica para ser seleccionado positivamente, a diferencia de otros supresores de tumores con funciones haploinsuficientes en el epitelio colónico^[37]. Sin embargo, nuestra inferencia de los sesgos de fijación se vio limitada por la escasez de mutaciones truncantes de APC detectadas, por lo que esta interpretación se beneficiaría de datos de mutación adicionales adecuados.

Las estimaciones previas sugieren que la tasa de transformación de una cripta que porta cualquier mutación conductora del cáncer es muy baja, aproximadamente una de cada tres millones, lo que dificulta la conciliación de la carga mutacional de los eventos de mutación conductora individuales con la incidencia del cáncer^[1]. Para los modelos evolutivos del CRC y los datos epidemiológicos, se necesitan entre tres y cinco mutaciones conductoras que ocurren simultáneamente para la iniciación^[6]–^[9]. Sin embargo, las criptas que portan dos o más mutaciones conductoras del CRC canónicas son extremadamente raras en el colon normal, lo que plantea preguntas sobre cómo surgen estas combinaciones durante la tumorigenesis temprana. Históricamente, la ubicuidad de la inactivación bialélica de APC en el CRC ha posicionado la pérdida de APC como el evento fundador en la mayoría de los tumores. En apoyo de esto, las simulaciones de la coocurrencia y el orden temporal de las mutaciones conductoras indican que si la selección para la pérdida bialélica de APC supera con creces los sesgos más modestos asociados con las mutaciones de TP53 o las mutaciones KRAS G12 monoalélicas, entonces los criptas nulos para APC tienen muchas más probabilidades de adquirir conductores adicionales a lo largo de la vida. Sin embargo, nuestras simulaciones también revelaron que en los individuos más jóvenes, donde la ventana de evolución clonal es más corta, surgen vías evolutivas alternativas, en las que un evento fundador KRAS es igualmente probable dentro de las criptas con triple mutación, lo que es particularmente relevante en el contexto del aumento de las tasas de incidencia del CRC de inicio temprano^[10].

Una de las principales incertidumbres al predecir el grado en que las vías evolutivas APC-primarias dominan la iniciación del CCRC es la verdadera magnitud y el momento del sesgo de selección conferido por la pérdida bialélica de APC. Las estimaciones actuales se basan en gran medida en modelos de ratón y en tejido de adenoma, y generalmente asumen que el comportamiento sesgado surge inmediatamente después de la pérdida del segundo alelo de APC^[20],^[25],^[36]. Sin embargo, los estudios en ratones sugieren que los fenotipos APC-nulos dependen de la densidad, de modo que se requiere un umbral crítico de criptas mutantes en APC para iniciar la formación de adenomas, y las bajas densidades no progresan y las altas densidades promueven la agregación y la iniciación del tumor^[38]. Estos hallazgos implican que la selección de clones APC-nulos puede depender del contexto en lugar de ser intrínseca. En apoyo de esto, un estudio reciente de secuenciación del genoma completo de criptas individuales de pacientes con poliposis adenomatosa familiar mostró que la inactivación somática del segundo alelo de APC casi siempre coincide con focos de criptas aberrantes (FCA) o pólipos adenomatosos y rara vez se encuentra en criptas morfológicamente normales^[39]. La detección morfológica consistente de la pérdida bialélica de APC sugiere que las distribuciones del tamaño de los FCA podrían utilizarse para cuantificar el sesgo de selección asociado a APC en el colon humano. Estas perspectivas cuantitativas serán esenciales para determinar si APC-primario es una vía obligatoria para el CCRC o una de las múltiples vías evolutivas que podrían depender de la edad.

Más allá de su paisaje genómico, existen otros factores fenotípicos y externos que determinan las limitaciones de la selección y el resultado neoplásico de los clones mutantes. Por ejemplo, la heterogeneidad en el resultado fenotípico de las mutaciones activadoras de KRAS documentadas aquí sugiere que existen diferencias estables e inherentes en las respuestas epiteliales a mutaciones idénticas, como su estado de linaje, lo que puede promover aún más su potencial transformador. Estos hallazgos reflejan los problemas de los esfuerzos de detección temprana que requieren la estratificación de las lesiones según su riesgo de transformación, para lo cual la evaluación patológica podría no ser suficiente. Se suele pensar que los pólipos hiperplásicos con mutaciones de KRAS y cambios morfológicos claros no progresan a CCRC, mientras que los clones mutantes de KRAS sin cambios morfológicos, como se observa aquí, presentan una firma transcriptómica que se asemeja a un subtipo específico de tumores de CCRC (CMS3/iCMS3). Esto podría proporcionar una base para probar estrategias de estratificación basadas en los fenotipos de KRAS mutantes, lo que podría facilitarse en un entorno clínico mediante la inmunodetección de REG4.

Métodos

Tejido humano

Se obtuvieron muestras de tejido de colon histológicamente normal de pacientes con cáncer (normal del cáncer) del Addenbrooke's Hospital Cambridge, Norwich University Hospital y St James University Hospital Leeds, con la aprobación completa del comité de ética de investigación local (REC 15/WA/0131 aprobado por el Comité de Ética de la Investigación de Gales 7, 06/Q0108/307 y 08/H0304/85, ambos aprobados por el Comité de Ética de la Investigación del Este de Inglaterra-Cambridge, así como 12/LO/1217 aprobado por el Comité de Ética de la Investigación de Londres-Bloomsbury, respectivamente), de acuerdo con la institución de la Autoridad de Investigación en Salud. Los especímenes de colectomía se fijaron en formaldehído tamponado neutro al 10%. De las áreas de tejido sin enfermedad macroscópicamente visible, se retiraron láminas de mucosa de los especímenes y se incluyeron en face en bloques de parafina (FFPE).

En total, el estudio contiene datos de 222 pacientes con una edad media de 68 años (23 a 91 años), 106 mujeres y 115 hombres. La cohorte total de pacientes se equilibró intencionalmente en cuanto al género y se incluyó un amplio rango de edades. Se ha obtenido el consentimiento del paciente para la publicación de identificadores indirectos, incluido el género y la edad.

Micro-seq

Las secciones se montaron en portaobjetos especiales con membrana de PET (Zeiss, 15511306) que se irradiaron con UV durante 30 minutos antes del corte. Los portaobjetos de disección micro láser (DML) se desparafinaron utilizando Leica Multistainer ST5020. La disección micro láser (DML) se realizó en el sistema de disección micro láser Leica LMD7000. Se seleccionaron y recogieron áreas de 200 criptas en tapas de tubos de 0,2 mm de radio y se resuspendieron en 30 μl de tampón de proteinasa K del kit de ADN Arcturus PicoPure (ThermoFisher, KIT0103). La lisis se realizó durante 5 horas a 65 °C, seguida de una incubación de 10 minutos a 95 °C. Después de la DML, se escanearon los portaobjetos y se anotó la ubicación de cada parche.

AceStrip- y subAceStrip-seq

Las secciones se montaron en láminas de acetato de 240 μm que se pre-recubrieron con 0,01% de poli-L-lisina (P8920-100ML) durante 5 minutos con agitación, seguido de un enjuague con agua destilada y se dejaron secar por completo. Las láminas de acetato se irradiaron con UV durante 30 minutos antes del corte. Después del corte, se almacenaron a 4 °C en un soporte personalizado.

La desparafinación se realizó manualmente utilizando un soporte personalizado (Figura suplementaria 6). Las secciones de acetato se sumergieron en xileno durante 10 minutos y luego en una serie de concentraciones decrecientes de etanol (100, 75 y 50%) durante 10 minutos cada una y, finalmente, se enjuagaron en agua destilada y se dejaron secar por completo. Las secciones sobre acetato se colocaron sobre una plataforma de plástico con una plantilla de metal encima. La plantilla se diseñó e imprimió en 3D para contener ranuras guía con centros separados por 2 mm. Las secciones se cortaron en tiras verticalmente, que a veces se dividieron en dos según la cantidad de epitelio en cada tira, con el objetivo de tener alrededor de 1000 criptas en cada tira para AceStrip-seq. Se tomaron imágenes de las tiras para anotar su ubicación y relación espacial. Para subAceStrip-seq, se utilizaron secciones seriales de AceStrip-seq, de modo que las tiras que contenían una mutación de interés se cortaron en trozos más pequeños, que contenían alrededor de 200 criptas cada uno.

Se utilizó el kit de ADN QIAmp FFPE Tissue (QIAGEN, 56404) según las instrucciones del fabricante para la extracción de ADN (eluido en 40 μl de agua).

Los datos del estudio incluyen datos de un total de 222 pacientes con una edad media de 68 años (23 a 91 años), 106 mujeres y 115 hombres. La cohorte total de pacientes se equilibró intencionalmente en cuanto al género y se incluyó un amplio rango de edades. Se ha obtenido el consentimiento del paciente para la publicación de identificadores indirectos, incluido el género y la edad.

Sección-seq

La secuenciación de ADN de sección completa se realizó utilizando espirales de ADN de secciones completas. Se utilizó el kit de ADN QIAmp FFPE Tissue (QIAGEN, 56404) según las instrucciones del fabricante para la extracción de ADN (eluido en 20 μl de agua).

Preparación de la biblioteca y secuenciación de última generación

Se diseñaron dos conjuntos de cebadores para abarcar un área de interés de un punto caliente en diferentes secuencias utilizando Primer3 para la amplificación del objetivo (Tabla suplementaria 3). Los adaptadores Illumina CS1 y CS2 se añadieron al principio de los cebadores directos o inversos, respectivamente. Los grupos de multiplexación de cebadores se enumeran en la Tabla suplementaria 4.

Los componentes de la reacción de PCR se describen en la Tabla suplementaria 5. La ciclación de PCR se realizó a 95 °C durante 2 minutos para un ciclo, seguido de 36 ciclos a 95 °C durante 10 segundos, 60 °C durante 10 segundos y 72 °C durante 15 segundos. El ciclo final fue seguido de una extensión de 5 minutos a 72 °C.

Las muestras se codificaron por barras utilizando el sistema de PCR Fast Start High Fidelity (Roche, Basilea, Suiza) según el protocolo del proveedor. Después de la agrupación y la purificación mediante el kit Clean & Concentrator (Zymo Research) y la selección del tamaño mediante PippinBlue (Sage Science), las muestras se secuenciaron utilizando una secuenciación apareada de 150 pares de bases (pb) con un 15% de PhiX en la plataforma Illumina NovaSeq (SP).

Llamada de mutaciones

La llamada de mutaciones se realizó utilizando RePlow^[40] utilizando la cobertura de amplicones de replicación dual. En resumen, RePlow explota las preparaciones de biblioteca de replicación para separar la contribución de los errores de fondo (en la etapa de preparación de la biblioteca) de los errores de secuenciación (en la etapa de secuenciación de la biblioteca) utilizando un modelo estadístico, que calcula la relación de probabilidad total para los candidatos a variantes en comparación con los errores de fondo en todas las replicaciones simultáneamente. Al hacerlo, se construyen los valores de frecuencia alélica ajustados (VAF), que se utilizaron como entrada para la inferencia clonal como se describe en la sección correspondiente de "Métodos", restando la contribución de los errores de secuenciación, mientras que los perfiles de error de fondo para el modelo de error se calculan de forma independiente para los seis tipos de sustitución de pares de bases (A > C, A > G, A > T, C > A, C > G, C > T) en las regiones objetivo. Para las llamadas con una relación de probabilidad positiva devuelta por RePlow, se aplicó un paso de filtrado estricto adicional, en el que solo se retuvieron las variantes que también tenían una relación de probabilidad positiva en ambas replicaciones individualmente. Esto fue particularmente importante para excluir las llamadas falsas positivas a bajos valores de VAF, que corresponden a un pequeño número de criptas mutantes sobre un fondo de tipo salvaje grande. Para subAceStrip-seq, las mutaciones GCC > T representaron 2/3 de todas las llamadas, destacando al realizar el análisis de la firma mutacional y siendo únicas para esa escala de secuenciación. Todas las llamadas en este contexto de trinucleótido se eliminaron ya que se consideraron artefactuales (posiblemente porque las secciones se dejaron sobre acetato durante más tiempo).

Algoritmo de teselación para la llamada de clones

Para asignar de forma eficiente múltiples llamadas de la misma mutación realizadas del mismo paciente a clones individuales, calculamos una teselación de Voronoi de las secciones de tejido utilizando las coordenadas espaciales de los sitios seleccionados para Micro-seq. De esta manera, cada tesela de Voronoi captura el área más grande de tejido no muestreado asignada a un sitio Micro-seq determinado. Múltiples llamadas de la misma mutación que se asignan a teselas de Voronoi adyacentes se fusionaron inicialmente en el mismo clon (mediante un algoritmo de búsqueda en profundidad), mientras que aquellas que estaban espacialmente separadas por teselas de Voronoi vacías (para esa mutación) se consideraron eventos clonales distintos. Los recuentos de los clones fusionados generados de esta manera se utilizaron como entrada para dN/dS. En comparación, se consideraron todas las posibles particiones de criptas mutantes entre clones, que son consistentes con los valores de VAF observados para estas llamadas fusionadas, durante la inferencia de la dinámica clonal, como se describe a continuación.

Firma de mutación

Se derivó una firma de mutación del exoma colónico humano mediante la recopilación de los 20 581 SNV identificados dentro de la región de codificación genómica de criptas colónicas humanas normales individuales por Lee-Six et al. Los contextos de trinucleótido de los SNV se transformaron al formato Sanger convencional y se sumaron los recuentos. Este proceso se repitió para las llamadas de mutación de Micro-seq. Luego, se normalizaron los recuentos de mutación de Micro-seq y Lee-Six et al. mediante la frecuencia de destino de trinucleótido dentro de sus respectivos paneles de secuenciación. La normalización proporciona una frecuencia relativa de cada SNV de contexto de trinucleótido que ocurre y permite una comparación directa entre nuestra firma de mutación y la derivada de Lee-Six et al., teniendo en cuenta las diferencias en la abundancia de tri-nucleótidos entre el panel de secuenciación relativamente pequeño utilizado aquí (1180 pb) y el exoma completo.

Cuantificación de la selección utilizando dN/dS

La última versión de la implementación de máxima verosimilitud de dN/dS descrita originalmente en la ref. ^[24] y disponible en el paquete R dNdScv (https://github.com/im3sanger/dndscv), se aplicó para identificar genes bajo selección positiva. Se construyó un objeto de referencia personalizado a partir del ensamblaje GRCh38/hg38 mediante la subselección de la matriz de consecuencias de sustitución dependiente del contexto de trinucleótido al conjunto de mutaciones posibles basadas en la región de 1180 pb del genoma seleccionada para la secuenciación dirigida en este estudio. Los recuentos de mutación utilizados para obtener los valores de dN/dS informados en el texto principal son de las llamadas de clon fusionado que utilizan una teselación de Voronoi como se describe en la sección anterior de "Métodos", pero repetimos el análisis con el supuesto más conservador de solo una mutación/clon única por paciente para obtener resultados comparables (coeficientes de dN/dS > 2 y q ajustado = 4,0 × 10−7 y 2,0 × 10−4 para las mutaciones de FBXW7 y KRAS con cambio de sentido, respectivamente). Las mutaciones de CTNNB1 también se incluyeron al ejecutar dN/dS, pero no se representan aquí (no se seleccionaron positivamente).

Cálculo de la relación estroma:epitelio

La relación cripta/estroma se calculó utilizando secciones seriadas de H&E y Halo (Indica Labs) para el análisis de imágenes (Figura suplementaria 6f, g). Se utilizó un clasificador cripta/estroma (bosque aleatorio) y un algoritmo de segmentación celular (IHC multiplexada) para calcular el número de células en el epitelio y el estroma. Se calculó una única relación estroma/epitelio por sección para los datos de Micro-seq y se calculó una relación estroma/epitelio individual por franja para Strip-seq. Para Section-seq, se asumió una relación estroma/epitelio de 1:1. Esto se utilizó para calcular el tamaño de parche inferido en función de la frecuencia alélica de la variante (VAF) de la siguiente manera: [2 × VAF × (1 + estroma/epitelio) × cripta N].

Inferencia de la dinámica clonal a partir de datos de secuenciación

Anteriormente, se obtuvieron estimaciones de las tasas de fijación y fisión conferidas por mutaciones específicas de KRAS (en las posiciones de aminoácidos G12 y G13) utilizando un modelo estadístico para la carga mutacional observada basado en datos de secuenciación dirigida ^[3], ^[4]. En ese modelo estadístico original, todas las criptas mutadas detectadas en una muestra (calculadas a partir de los valores de VAF observados) se atribuyen a un solo clon. Sin embargo, en la práctica, cualquier combinación de clones individuales con un número total de criptas mutadas igual al observado puede explicar los datos de manera equivalente, una posibilidad que consideramos en el modelo estadístico actualizado que se describe aquí.

A partir de la distribución del tamaño de los clones ^[4], la probabilidad de un clon mutante con n criptas mutantes a la edad t es donde μ es la tasa de conversión monoclonal (tasa de fijación) de las criptas mutantes y ν es la tasa de fisión asociada a esa mutación. La probabilidad de que una cripta no esté mutada es

Dado un número observado de criptas mutantes n en una muestra con k clones totales, la probabilidad de una partición particular de n en k clones mutantes (con n_i el número de criptas mutantes en el clon i) es donde N es el número de clones en la partición con n_i y k con i elegidos para normalizar las probabilidades multinominales para que sumen uno. Tenga en cuenta que la partición asegura que Σn_i = n y Σk = k para todos los i si solo hay un clon con n criptas mutantes, por lo que la suposición original utilizada en ^[4] puede considerarse como un caso especial de este modelo más general.

Para calcular la probabilidad total de observar una carga mutacional de n criptas mutantes, se deben sumar las probabilidades de todas las posibles particiones de n. La enumeración de particiones enteras conduce a una explosión combinatoria a medida que aumenta n, pero sumar sobre particiones durante la inferencia de parámetros es factible desde el punto de vista computacional, siempre que n permanezca suficientemente pequeño. En nuestros datos de Micro-seq, todos los números de criptas mutadas observadas, excepto uno (una mutación KRAS G12C encontrada en 15 criptas), fueron menores o iguales a 15, que, por lo tanto, se utilizó como tamaño máximo de clon durante la inferencia de la dinámica clonal. Como se describe para el método de inferencia original en Nicholson et al., se calculó un número crítico de criptas mutantes para las muestras en las que no se detectaron mutaciones para tener en cuenta la probabilidad de que el número total de criptas mutantes esté por debajo del límite de detección de la secuenciación. Los valores previos utilizados aquí para μ y ν son los mismos que en Nicholson et al.

Dado que el modelo estadístico actual tiene en cuenta, pero no puede resolver, la ambigüedad con respecto a la partición precisa que da lugar al número observado de criptas mutantes en cada muestra, los valores inferidos para μ y ν deben interpretarse como efectivos, ya que son los más consistentes con los datos observados y sus intervalos de credibilidad del 95% abarcan todas las posibilidades. Las simulaciones confirmaron que el modelo basado en particiones genera estimaciones razonables para las tasas de fisión y fijación efectivas en escenarios extremos donde las tasas de fisión reales son comparativamente bajas y las tasas de fijación reales son muy altas, lo que el modelo original no logra hacer al asignar todas las criptas mutantes a un solo clon (Figura suplementaria 7). Al aplicar el modelo a datos reales de grupos de variantes que se supone que comparten tasas de fisión comparables, consideramos dos regímenes diferentes: ya sea que al grupo se le asignara un único valor compartido de ν o que a cada variante se le asignara su propio ν de forma independiente, como se describe en el texto principal. ν se utilizó para generar figuras con distribuciones ilustrativas del tamaño de los clones para diferentes grupos de mutaciones utilizando los valores inferidos de μ y ν con una edad promedio del paciente de 60 años.

Simulaciones de carga mutacional

Las tasas de fijación y fisión inferidas se utilizaron para calcular la carga mutacional estimada y trazarla para un grupo de mutaciones en el colon normal a una edad determinada, como se describe en la referencia ^[37] según la fórmula:

Aquí, P es la distribución de la probabilidad de un parche de criptas mutadas de tamaño n a la edad t, como se describe en la sección anterior de "Métodos". La carga mutacional de una coocurrencia triple se predijo entonces asumiendo la independencia de los tres eventos individuales multiplicando las curvas de carga mutacional correspondientes. Para expresar la carga mutacional en términos del número de criptas mutantes por individuo, utilizamos N = 42 millones, asumiendo un número total de 42 millones de criptas en el colon humano promedio.

Para simular las cargas mutacionales en ausencia de selección, asignamos la tasa de fisión inferida para las mutaciones sinónimas a todos los grupos de mutaciones y utilizamos la estimación para las mutaciones sinónimas en el contexto del trinucleótido GCG>T para FBXW7 R465C. Las observaciones experimentales de la carga mutacional obtenidas de AceStrip-seq y Section-seq se calcularon para cada variante dividiendo el número estimado de criptas mutantes por el número total de criptas examinadas a una edad determinada, excluyendo a los pacientes menores de 50 años para evitar puntos de datos anómalos asociados con una carga mutacional excepcionalmente alta a una edad temprana. Los errores cuadráticos medios se calcularon en función de la distancia euclidiana media entre las cargas mutacionales observadas y predichas en ausencia o presencia de selección, siendo la relación entre estas la medida del rendimiento (la relación del error cuadrático medio (RRMSE) = 1,0 si el modelo con selección describe los datos tan bien como el modelo sin selección, mientras que RRMSE < 1,0 cuando el modelo con selección funciona mejor).

Simulaciones y comparaciones con las tasas de incidencia de cáncer

El modelo matemático de Paterson et al. se modificó para incorporar las tasas de fijación y fisión de los clones inferidas utilizando los grupos de mutaciones truncantes de APC, KRAS G12 y las mutaciones impulsoras seleccionadas de TP53 detectadas en nuestro estudio. En resumen, el modelo original se parametrizó utilizando las tasas de mutación de novo específicas del gen estimadas en función del número de posiciones que generan mutaciones inactivadoras (para APC y TP53) y activadoras (para KRAS) en el gen correspondiente, combinadas con datos experimentales de una mezcla de estudios en ratones y humanos. Aquí, las tasas de fijación y fisión de las mutaciones monoalélicas en estos genes se tomaron directamente de nuestra dinámica clonal inferida, y el modelo modificado debe interpretarse, por lo tanto, como uno que restringe aún más el conjunto de posibles mutaciones a aquellas detectadas mediante Micro-seq. Los parámetros no medidos correspondientes a la tasa de pérdida de heterocigosidad mantuvieron los mismos valores utilizados en el modelo original, mientras que las tasas de fijación y fisión de los clones con inactivación bialélica de APC (y, por lo tanto, criptas con inactivación bialélica de APC y activación de KRAS) se consideraron ya sea: idénticas a las de la inactivación monoalélica de APC (modelo sin sesgo de APC); o sin cambios con respecto al modelo original (modelo con sesgo de APC). Ambos modelos se simularon utilizando la configuración predeterminada y se optimizaron para el procesamiento en paralelo en 120 núcleos de HPC siguiendo la documentación del repositorio de código (https://github.com/chaypaterson/CRC_Initiation), primero utilizando la versión tau-leaping del algoritmo de simulación estocástica para calcular las probabilidades de criptas mutantes a una edad determinada y luego en modo de linaje para rastrear el orden de los diferentes eventos de mutación, como se describe en la publicación original ^[25].

Las probabilidades de criptas mutantes a una edad determinada se compararon con los datos de los siguientes estudios de cbioportal (https://www.cbioportal.org/)^[41] que habían cubierto las regiones incluidas en el panel de amplicones que comprenden el conjunto de datos de cáncer colorrectal: cáncer de colon (Sidra-LUMC AC-ICAM, Nat Med 2023), adenocarcinoma colorrectal (DFCI, Cell Reports 2016), adenocarcinoma colorrectal (Genentech, Nature 2012), adenocarcinoma colorrectal (TCGA, Firehose Legacy), cáncer colorrectal (CAS Shanghai, Cancer Cell 2020), cáncer colorrectal (MSK, Gastroenterology 2020), "cáncer colorrectal (MSK, JNCI 2021), cáncer colorrectal metastásico (MSK, Cancer Cell 2018), cáncer rectal (MSK, Nature Medicine 2022), cáncer rectal (MSK, Nature Medicine 2019). La frecuencia de una coocurrencia triple en el CRC se calculó como la fracción de muestras que contienen al menos una mutación de cada uno de los tres conjuntos de mutaciones correspondientes, dividida por el número total de muestras (N = 6075), con intervalos de confianza del 95% obtenidos utilizando el método de Clopper-Pearson. Esta frecuencia se multiplicó luego por el riesgo de por vida de CRC (5%) para estimar la incidencia de por vida de CRC con esa coocurrencia triple en la población.

Inmunohistoquímica/inmunofluorescencia

La IHC se realizó como se describió anteriormente ^[4]. En resumen, las secciones de parafina se desparafinaron y se rehidrataron utilizando lavados secuenciales de xileno y alcoholes graduados, seguidos de la recuperación del antígeno mediada por calor en un tampón de citrato de trisodio de 10 mM (pH 6,0). Las secciones se incubaron en peróxido de hidrógeno al 3% en metanol durante 15 minutos y se bloquearon con suero de burro al 10% durante 30 minutos. Luego, las secciones se incubaron durante la noche a 4 °C con anticuerpos primarios contra REG4 (HPA, HPA046555, 1:1000, policlonal) y MUC5AC (Abcam, ab198294, 1:200, clon EPR16904). Después de lavar las láminas, se incubaron durante 40 minutos a temperatura ambiente con anticuerpo secundario biotin-SP-conjugado AffiniPure de burro anti-conejo (Jackson ImmunoResearch; 1:500 en PBS-T), seguido de la incubación con el reactivo VECTASTAIN® Elite® ABC (Vector Laboratories) durante 40 minutos. La detección de la señal se realizó utilizando un sistema de cromógeno líquido DAB+ (Dako). Las secciones se contracolorearon con hematoxilina, se deshidrataron, se aclararon y se montaron.

El protocolo de IHC se adaptó para inmunofluorescencia sin los pasos de peróxido de hidrógeno y DAB y, en su lugar, utilizando anticuerpos secundarios de fluoróforo Alexa conjugados (Thermo Fisher, 1:200) y DAPI (1 μg/mL). Se utilizaron los anticuerpos primarios MUC5AC (Abcam, ab198294, 1:100, clon EPR16904) y CDX2 (Sigma, AMAB91828, 1:500, clon CL12974). Las láminas se escanearon en el sistema de imagenología cuantitativa automatizada Vectra Polaris de Akoya Biosciences (Perkin Elmer).

El anticuerpo REG4 ha sido validado por Human Protein Atlas. El anticuerpo MUC5AC ha sido validado por Abcam y ha sido citado en más de 10 publicaciones. El anticuerpo CDX2 ha sido validado por Sigma utilizando RNAseq ortogonal.

Acumulación de edad e inferencia basada en IHC de REG4

Se utilizó el paquete RHClones (https://github.com/ElEd2/RHClones) para estimar la pendiente relacionada con la edad de los clones detectados (N = 291). Se utilizó el mismo paquete para inferir la dinámica clonal basada en los datos de IHC de REG4.

LCM para muestras positivas para REG4

Se cortaron secciones seriadas de muestras con regiones positivas para REG4 en secciones LCM, como se describió anteriormente (sección de Micro-seq). Las regiones positivas para REG4, junto con algunas criptas WT adyacentes para aumentar la entrada de material genómico, se cortaron y se secuenciaron utilizando la secuenciación de próxima generación para las mutaciones de KRAS (G12/G13, A146, Q61) utilizando los cebadores mencionados anteriormente (Tabla suplementaria 3). Se incluyeron muestras WT/negativas para REG4 para cada paciente para excluir cualquier posible SNP. Para la llamada de mutaciones, se utilizó un ensamblador de pares finales para secuencias de Illumina, PANDAseq, para fusionar los pares de lecturas correspondientes en un amplicon in silico ^[42]. Se identificaron las diez especies diferentes principales que comenzaron y terminaron con la secuencia de cebador correcta y se calculó su frecuencia. Las muestras mutadas contenían especies WT junto con especies mutantes de KRAS.

Procesamiento de la lámina Visium y preparación de la biblioteca

Se cortaron secciones seriadas con los clones mutantes de interés identificados a partir de MicroSeq y AceStripseq, colocando la ROI (región de interés) en el centro de la lámina de vidrio. El procesamiento de las láminas y la preparación de las bibliotecas se realizaron según el protocolo 10X Visium Cytassist FFPE.

Análisis de datos de Visium

La alineación de datos se realizó utilizando Spaceranger (10X). Se utilizó LoupeBrowser para dibujar alrededor de las regiones mutantes de KRAS y las regiones WT (tipo salvaje) correspondientes, y para extraer los códigos de barras. Los puntos con UM1 < 1000 genes < 900 y porcentaje de genes mitocondriales > 20% se eliminaron del análisis. Se realizó un análisis de componentes principales, un agrupamiento no supervisado, un análisis de expresión génica diferencial y una reducción de dimensionalidad utilizando la biblioteca de experimentos espaciales en R.

La firma mutante de KRAS se generó compilando los 400 genes con mayor expresión diferencial de cada uno de los 14 clones y seleccionando los genes comunes a al menos dos clones (N = 628 genes únicos). Se realizó GSEA utilizando las bibliotecas cluster Profiler y msigdbr en R.

Dado que identificamos que R465C se asocia con un fuerte sesgo de selección, utilizamos la información de localización espacial de los clones proporcionada por Micro-seq y AceStrip-seq para realizar un análisis de reducción de dimensionalidad en las regiones más pequeñas que contienen los clones mutantes de FBXW7 R465C. No apareció ningún grupo distinto en estas regiones, lo que probablemente refleja la falta de alteraciones transcripcionales asociadas con la pérdida de un alelo supresor tumoral y/o el pequeño tamaño de los clones, lo que impidió el agrupamiento.

Pruebas estadísticas

Todas las pruebas estadísticas se realizaron utilizando las funciones predeterminadas disponibles en el lenguaje de programación R. Se realizó la prueba exacta de Fisher utilizando la función fisher.exact() con los valores de los parámetros predeterminados.

Resumen del informe

Se puede encontrar información adicional sobre el diseño de la investigación en el Resumen del informe de Nature Portfolio vinculado a este artículo.

Tejido humano

Se obtuvieron muestras de tejido de colon histológicamente normal de pacientes con cáncer (tejido normal del cáncer) del Addenbrooke’s Hospital Cambridge, Norwich University Hospital y St James University Hospital Leeds, bajo la aprobación completa del comité de ética de investigación local (REC 15/WA/0131 aprobado por el Comité de Ética de la Investigación de Gales 7, 06/Q0108/307 y 08/H0304/85, ambos aprobados por el Comité de Ética de la Investigación del Este de Inglaterra-Cambridge Este, así como 12/LO/1217 aprobado por el Comité de Ética de la Investigación de Londres-Bloomsbury, respectivamente), de acuerdo con la Health Research Authority. Los especímenes de colectomía se fijaron en formaldehído tamponado neutro al 10%. De las áreas de tejido sin enfermedad macroscópicamente visible, se retiraron las láminas mucosas de los especímenes y se incluyeron en face en bloques de parafina (FFPE).

En total, el estudio contiene datos de 222 pacientes con una edad promedio de 68 años (23-91 años), 106 mujeres y 115 hombres. La cohorte total de pacientes se equilibró intencionalmente en cuanto al género y se incluyó un amplio rango de edades. Se obtuvo el consentimiento del paciente para la publicación de identificadores indirectos, incluido el género y la edad.

Micro-seq

Las secciones se montaron en láminas especiales con membrana de PET (Zeiss, 15511306) que se irradiaron con UV durante 30 minutos antes del corte. Las láminas de LCM (microdisección láser de captura) se desparafinaron utilizando Leica Multistainer ST5020. La microdisección láser de captura (LCM) se realizó en el sistema de microdisección láser Leica LMD7000. Se seleccionaron y recolectaron áreas de 200 criptas en tapas de tubos de 0,2 mm de radio y se resuspendieron en 30 μl de tampón de Proteinasa K del kit de ADN Arcturus PicoPure (ThermoFisher, KIT0103). La lisis se realizó durante 5 horas a 65 °C, seguida de una incubación de 10 minutos a 95 °C. Después de la LCM, se escanearon las láminas y se anotó la ubicación de cada parche.

AceStrip- y subAceStrip-seq

Las secciones se montaron en láminas de acetato de 240 μm que se habían pre-recubierto con 0,01% de poli-L-lisina (P8920-100ML) durante 5 minutos con agitación, seguido de un enjuague con agua destilada y se dejaron secar por completo. Las láminas de acetato se irradiaron con UV durante 30 minutos antes del corte. Después del corte, se almacenaron a 4 °C en un soporte personalizado.

La desparafinación se realizó manualmente utilizando un soporte personalizado (Figura suplementaria 6). Las secciones de acetato se sumergieron en xileno durante 10 minutos y luego en una serie de concentraciones decrecientes de etanol (100, 75 y 50%) durante 10 minutos cada una, y finalmente se enjuagaron en agua destilada y se dejaron secar por completo. Las secciones sobre acetato se colocaron sobre una plataforma de plástico con una plantilla de metal en la parte superior. La plantilla se diseñó e imprimió en 3D para contener ranuras guía con centros separados por 2 mm. Las secciones se cortaron en tiras verticalmente, que a veces se dividieron en dos según la cantidad de epitelio en cada tira, con el objetivo de tener alrededor de 1000 criptas en cada tira para AceStrip-seq. Se tomaron imágenes de las tiras para registrar su ubicación y relación espacial. Para subAceStrip-seq, se utilizaron secciones seriadas de AceStrip-seq, donde las tiras que contenían una mutación de interés se cortaron en trozos más pequeños, que contenían alrededor de 200 criptas cada uno.

Se utilizó el kit de ADN QIAmp FFPE Tissue (QIAGEN, 56404) según las instrucciones del fabricante para la extracción de ADN (eluido en 40 μl de agua).

Section-seq

La secuenciación de ADN de sección completa se realizó utilizando espirales de ADN de secciones completas. Se utilizó el kit de ADN QIAmp FFPE Tissue (QIAGEN, 56404) según las instrucciones del fabricante para la extracción de ADN (eluido en 20 μl de agua).

Preparación de bibliotecas y secuenciación de última generación

Se diseñaron dos conjuntos de cebadores para abarcar un área de interés en un punto caliente en diferentes secuencias utilizando Primer3 para la amplificación del objetivo (Tabla suplementaria 3). Los adaptadores CS1 y CS2 de Illumina se agregaron al principio de los cebadores directos o inversos, respectivamente. Los grupos de multiplexación de cebadores se enumeran en la Tabla suplementaria 4.

Los componentes de la reacción de PCR se describen en la Tabla suplementaria 5. La ciclación de PCR se realizó a 95 °C durante 2 minutos para un ciclo, seguido de 36 ciclos a 95 °C durante 10 segundos, 60 °C durante 10 segundos y 72 °C durante 15 segundos. El ciclo final fue seguido por una extensión de 5 minutos a 72 °C.

Las muestras se codificaron por barras utilizando el sistema de PCR Fast Start High Fidelity (Roche, Basilea, Suiza) según el protocolo del proveedor. Después de la agrupación y la purificación mediante el kit Clean & Concentrator (Zymo Research) y la selección del tamaño mediante PippinBlue (Sage Science), las muestras se secuenciaron utilizando secuenciación apareada de 150 pares de bases (bp) con 15% de PhiX en el sistema Illumina NovaSeq (SP).

Llamada de mutaciones

La llamada de mutaciones se realizó utilizando RePlow^[40] utilizando la cobertura de amplicones de replicación dual. En resumen, RePlow explota las preparaciones de bibliotecas de replicación para separar la contribución de los errores de fondo (en la etapa de preparación de la biblioteca) de los errores de secuenciación (en la etapa de secuenciación de la biblioteca) utilizando un modelo estadístico, que calcula la relación de probabilidad total para los candidatos a variantes en comparación con los errores de fondo en todas las replicaciones simultáneamente. Al hacerlo, se construyen VAF ajustados (que se utilizaron como entrada para la inferencia clonal como se describe en la sección correspondiente de 'Métodos') restando la contribución de los errores de secuenciación, mientras que los perfiles de error de fondo para el modelo de error se calculan de forma independiente para los seis tipos de sustitución de pares de bases (A > C, A > G, A > T, C > A, C > G, C > T) en las regiones objetivo. Para las llamadas con una relación de probabilidad positiva devuelta por RePlow, se aplicó un paso de filtrado estricto adicional, en el que solo se retuvieron las variantes que también tenían una relación de probabilidad positiva en ambas replicaciones individualmente. Esto fue particularmente importante para excluir las llamadas falsas positivas a valores bajos de VAF, lo que corresponde a un pequeño número de criptas mutantes sobre un gran fondo de tipo salvaje. Para subAceStrip-seq, las mutaciones GCC > T representaron 2/3 de todas las llamadas, destacando al realizar el análisis de la firma mutacional y siendo únicas de esa escala de secuenciación. Todas las llamadas en este contexto de trinucleótido se eliminaron ya que se consideraron artefactuales (posiblemente porque las secciones se dejaron sobre acetato durante más tiempo).

Algoritmo de teselación para la llamada de clones

Para asignar de forma parsimoniosa múltiples llamadas de la misma mutación realizadas del mismo paciente a clones individuales, calculamos una teselación de Voronoi de las secciones de tejido utilizando las coordenadas espaciales de los sitios seleccionados para Micro-seq. De esta manera, cada tesela de Voronoi captura el área más grande de tejido no muestreado asignada a un sitio Micro-seq determinado. Múltiples llamadas de la misma mutación que se mapean a teselas de Voronoi adyacentes se fusionaron inicialmente en el mismo clon (a través de un algoritmo de búsqueda en profundidad), mientras que aquellas que estaban espacialmente separadas por teselas de Voronoi vacías (para esa mutación) se consideraron eventos clonales distintos. Los recuentos de los clones fusionados generados de esta manera se utilizaron como entrada para dN/dS. En comparación, se consideraron todas las posibles particiones de criptas mutantes entre clones, coherentes con los valores de VAF observados para estas llamadas fusionadas, durante la inferencia de la dinámica clonal, como se describe a continuación.

Firma de mutación

Se derivó una firma de mutación del exoma colónico humano mediante la recopilación de los 20.581 SNV identificados dentro de la región de codificación genómica de criptas colónicas humanas normales individuales por Lee-Six et al. Los contextos de trinucleótido de los SNV se transformaron al formato Sanger convencional y se sumaron los recuentos. Este proceso se repitió para las llamadas de mutación de Micro-seq. Luego, se normalizaron los recuentos de mutación de Micro-seq y los SNV de Lee-Six et al. mediante la frecuencia de destino de trinucleótido dentro de sus respectivos paneles de secuenciación. La normalización proporciona una frecuencia relativa de cada SNV de contexto de trinucleótido que ocurre y permite una comparación directa entre nuestra firma de mutación y la derivada de Lee-Six et al., teniendo en cuenta las diferencias en la abundancia de tri-nucleótidos entre el panel de secuenciación relativamente pequeño utilizado aquí (1180 pb) y el exoma completo.

Cuantificación de la selección utilizando dN/dS

La última versión de la implementación de máxima verosimilitud de dN/dS descrita originalmente en la ref. ^[24] y disponible en el paquete R dNdScv (https://github.com/im3sanger/dndscv), se aplicó para identificar genes bajo selección positiva. Se construyó un objeto de referencia personalizado a partir del ensamblaje GRCh38/hg38 mediante la subselección de la matriz de consecuencias de sustitución dependiente del contexto de trinucleótido al conjunto de posibles mutaciones en función de la región de 1180 pb del genoma seleccionada para la secuenciación dirigida en este estudio. Los recuentos de mutaciones utilizados para obtener los valores de dN/dS informados en el texto principal provienen de las llamadas de clones fusionados que utilizan una teselación de Voronoi como se describe en la sección anterior de 'Métodos', pero repetimos el análisis con la suposición más conservadora de solo una mutación/clon única por paciente para obtener resultados comparables (coeficientes de dN/dS > 2 y q ajustado = 4,0 × 10−7 y 2,0 × 10−4 para las mutaciones de FBXW7 y KRAS, respectivamente). Las mutaciones de CTNNB1 también se incluyeron al ejecutar dN/dS, pero no se representan aquí (no se seleccionaron positivamente).

Cálculo de la relación estroma:epitelio

La relación cripta a estroma se calculó utilizando secciones seriadas de H&E y Halo (Indica Labs) para el análisis de imágenes (Figura suplementaria 6f, g). Se utilizó un clasificador de criptas a estroma (bosque aleatorio) y un algoritmo de segmentación celular (Multiplex IHC) para calcular el número de células en el epitelio y el estroma. Se calculó una única relación estroma:epitelio por sección para los datos de Micro-seq y se calculó una relación estroma:epitelio individual por tira para Strip-seq. Para Section-seq, asumimos una relación estroma:epitelio de 1:1. Esto se utilizó para calcular el tamaño del parche inferido en función de la frecuencia alélica de la variante VAF de la siguiente manera: [2 × VAF × (1+ estroma/epitelio) × cripta N].

Inferencia de la dinámica clonal a partir de los datos de secuenciación

Anteriormente, las estimaciones de las tasas de fijación y fisión conferidas por mutaciones específicas de KRAS (en las posiciones de aminoácidos G12 y G13) se obtuvieron utilizando un modelo estadístico para la carga mutacional observada, basado en datos de secuenciación dirigida^[3],^[4]. En ese modelo estadístico original, todos los criptas mutados detectados en una muestra (calculados a partir de los valores de VAF observados) se atribuyen a un solo clon. Sin embargo, en la práctica, cualquier combinación de clones individuales con un número total de criptas mutadas igual al observado puede explicar los datos de manera equivalente, una posibilidad que se tiene en cuenta en el modelo estadístico actualizado que se describe aquí.

A partir de la distribución del tamaño de los clones^[4], la probabilidad de un clon mutado con n criptas mutadas a la edad t es donde r es la tasa de conversión monoclonal (tasa de fijación) de las criptas mutadas y f es la tasa de fisión asociada con esa mutación. La probabilidad de que una cripta no esté mutada es

Dado un número observado de criptas mutadas n en una muestra con k clones totales, la probabilidad de una partición particular de n en k clones mutados (con nᵢ el número de criptas mutadas en el clon i) es donde p es el número de clones en la partición con nᵢ y q es el número de clones en la partición con nᵢ, elegidos para normalizar las probabilidades multinomiales de modo que sumen uno. Tenga en cuenta que la partición asegura que nᵢ = n y k = k para todos los i si solo hay un clon con n criptas mutadas, por lo que la suposición original utilizada en ^[3] puede considerarse como un caso especial de este modelo más general.

Para calcular la probabilidad total de observar una carga mutacional de n criptas mutadas, se debe sumar las probabilidades de todas las particiones posibles de n. La enumeración de particiones enteras conduce a una explosión combinatoria a medida que n aumenta, pero la suma sobre particiones durante la inferencia de parámetros es factible desde el punto de vista computacional, siempre que n permanezca suficientemente pequeño. En nuestros datos de Micro-seq, todos los números observados de criptas mutadas, excepto uno (una mutación KRAS G12C encontrada en n criptas), fueron menores o iguales a n, que, por lo tanto, se utilizó como tamaño máximo de clon durante la inferencia de la dinámica de los clones. Como se describe para el método de inferencia original en Nicholson et al., se calculó un número crítico de criptas mutadas para las muestras en las que no se detectaron mutaciones, para tener en cuenta la probabilidad de que el número total de criptas mutadas esté por debajo del límite de detección de la secuenciación. Los valores previos utilizados aquí para r y f son los mismos que en Nicholson et al.

Dado que el modelo estadístico actual tiene en cuenta, pero no puede resolver, la ambigüedad con respecto a la partición precisa que da lugar al número observado de criptas mutadas en cada muestra, los valores inferidos para r y f deben interpretarse como efectivos, ya que son los más consistentes con los datos observados y sus intervalos de credibilidad del 95% abarcan todas las posibilidades. Las simulaciones confirmaron que el modelo basado en particiones genera estimaciones razonables para las tasas de fisión y fijación efectivas en escenarios extremos donde las tasas de fisión reales son comparativamente bajas y las tasas de fijación reales son muy altas, lo que el modelo original no logra hacer al asignar todas las criptas mutadas a un solo clon (Figura suplementaria 7). Al aplicar el modelo a datos reales de grupos de variantes que se supone que comparten tasas de fisión comparables, se consideraron dos regímenes diferentes: ya sea que al grupo se le asignara un único valor compartido de f o que a cada variante se le asignara su propio f de forma independiente, como se describe en el texto principal. Se utilizó f para generar figuras con distribuciones ilustrativas del tamaño de los clones para diferentes grupos de mutaciones, utilizando los valores inferidos de r y f con una edad promedio del paciente de 60 años.

Simulaciones de la carga mutacional

Se utilizaron las tasas de fijación y fisión inferidas para calcular la carga mutacional estimada y trazarla para un grupo de mutaciones en el colon normal a una edad t dada, como se describe en la referencia ^[37] según la fórmula Aquí, P(n, t) es la distribución de la probabilidad de un parche de criptas mutadas de tamaño n a la edad t, como se describe en la sección anterior de "Métodos". La carga mutacional de una triple coocurrencia se predijo entonces asumiendo la independencia de los tres eventos individuales, multiplicando las curvas de carga mutacional correspondientes. Para expresar la carga mutacional en términos del número de criptas mutadas por individuo, se utilizó N = 42 millones, asumiendo un número total de 42 millones de criptas en el colon humano promedio.

Para simular las cargas mutacionales en ausencia de selección, se asignó la tasa de fisión inferida para las mutaciones sinónimas a todos los grupos de mutaciones y se utilizó la estimación de f para las mutaciones sinónimas en el contexto del trinucleótido GCG>T para FBXW7 R465C. Las observaciones experimentales de la carga mutacional obtenidas de AceStrip-seq y Section-seq se calcularon para cada variante dividiendo el número estimado de criptas mutadas por el número total de criptas analizadas a una edad determinada, excluyendo a los pacientes menores de 50 años para evitar puntos de datos anómalos asociados con una carga mutacional excepcionalmente alta a una edad temprana. Los errores cuadráticos medios se calcularon en función de la distancia euclidiana media entre las cargas mutacionales observadas y predichas en ausencia o presencia de selección, siendo la relación entre ambas una medida del rendimiento (la relación del error cuadrático medio (RRMSE) = 1,0 si el modelo con selección describe los datos tan bien como el modelo sin selección, mientras que RRMSE < 1,0 cuando el modelo con selección funciona mejor).

Simulaciones y comparaciones con las tasas de incidencia de cáncer

El modelo matemático de Paterson et al. se modificó para incorporar las tasas de fijación y fisión de los clones inferidas utilizando los grupos de APC truncantes, KRAS G12 y las mutaciones impulsoras seleccionadas de TP53 detectadas en nuestro estudio. En resumen, el modelo original se parametrizó utilizando las tasas de mutación de novo específicas del gen estimadas en función del número de posiciones que generan mutaciones inactivadoras (para APC y TP53) y activadoras (para KRAS) en el gen correspondiente, combinadas con datos experimentales de una mezcla de estudios en ratones y humanos. Aquí, las tasas de fijación y fisión para las mutaciones monoalélicas en estos genes se tomaron directamente de nuestra dinámica de clones inferida, y el modelo modificado debe interpretarse, por lo tanto, como uno que restringe aún más el conjunto de mutaciones posibles a aquellas detectadas mediante Micro-seq. Los parámetros no medidos correspondientes a la tasa de pérdida de heterocigosidad conservaron los mismos valores utilizados en el modelo original, mientras que las tasas de fijación y fisión de los clones con inactivación bialélica de APC (y, por lo tanto, criptas con inactivación bialélica de APC y activación de KRAS) se consideraron ya sea: idénticas a las de la inactivación monoalélica de APC (modelo sin sesgo de APC); o sin cambios con respecto al modelo original (modelo con sesgo de APC). Ambos modelos se simularon utilizando la configuración predeterminada y se optimizaron para el procesamiento en 120 núcleos de HPC, siguiendo la documentación del repositorio de código (https://github.com/chaypaterson/CRC_Initiation), primero utilizando la versión de salto tau del algoritmo de simulación estocástica para calcular las probabilidades de criptas mutadas en una edad determinada y luego en modo de linaje para rastrear el orden de los diferentes eventos de mutación, como se describe en la publicación original^[25].

Se compararon las probabilidades de criptas mutadas en una edad determinada con los datos de los siguientes estudios de cbioportal (https://www.cbioportal.org/)^[41] que habían cubierto las regiones incluidas en el panel de amplicones que comprenden el conjunto de datos de cáncer colorrectal: cáncer de colon (Sidra-LUMC AC-ICAM, Nat Med 2023), adenocarcinoma colorrectal (DFCI, Cell Reports 2016), adenocarcinoma colorrectal (Genentech, Nature 2012), adenocarcinoma colorrectal (TCGA, Firehose Legacy), cáncer colorrectal (CAS Shanghai, Cancer Cell 2020), cáncer colorrectal (MSK, Gastroenterology 2020), "cáncer colorrectal (MSK, JNCI 2021), cáncer colorrectal metastásico (MSK, Cancer Cell 2018), cáncer rectal (MSK, Nature Medicine 2022), cáncer rectal (MSK, Nature Medicine 2019). La frecuencia de una triple coocurrencia en el CRC se calculó como la fracción de muestras que contenían al menos una mutación de cada uno de los tres conjuntos de mutaciones correspondientes, dividida por el número total de muestras (N = 6075), con intervalos de confianza del 95% obtenidos utilizando el método de Clopper-Pearson. Esta frecuencia se multiplicó luego por el riesgo de por vida de CRC (5%) para estimar la incidencia de por vida de CRC con esa triple coocurrencia en la población.

Inmunohistoquímica/inmunofluorescencia

La IHC se realizó como se describió anteriormente^[4]. En resumen, las secciones de parafina se desparafinaron y se rehidrataron utilizando lavados secuenciales de xileno y alcoholes graduados, seguidos de la recuperación del antígeno mediada por calor en un tampón de citrato de trisodio de 10 mM (pH 6,0). Las secciones se incubaron en peróxido de hidrógeno al 3% en metanol durante 15 minutos y se bloquearon con suero de burro al 10% durante 30 minutos. Luego, las secciones se incubaron durante la noche a 4 °C con anticuerpos primarios contra REG4 (HPA, HPA046555, 1:1000, policlonal) y MUC5AC (Abcam, ab198294, 1:200, clon EPR16904). Después del lavado, las diapositivas se incubaron durante 40 minutos a temperatura ambiente con anticuerpo secundario AffiniPure conjugado con biotina-SP de burro anti-conejo (Jackson ImmunoResearch; 1:500 en PBS-T), seguido de la incubación con el reactivo VECTASTAIN® Elite® ABC (Vector Laboratories) durante 40 minutos. La detección de la señal se realizó utilizando un sistema de cromógeno líquido DAB+ (Dako). Las secciones se contracolorearon con hematoxilina, se deshidrataron, se aclararon y se montaron.

El protocolo de IHC se adaptó para inmunofluorescencia sin los pasos de peróxido de hidrógeno y DAB, y en su lugar se utilizaron anticuerpos secundarios conjugados con fluoróforo Alexa (Thermo Fisher, 1:200) y DAPI (1 μg/mL). Se utilizaron los anticuerpos primarios MUC5AC (Abcam, ab198294, 1:100, clon EPR16904) y CDX2 (Sigma, AMAB91828, 1:500, clon CL12974). Las diapositivas se escanearon en el sistema de imagenología cuantitativa automatizada Vectra Polaris de Akoya Biosciences (Perkin Elmer).

El anticuerpo REG4 ha sido validado por Human Protein Atlas. El anticuerpo MUC5AC ha sido validado por Abcam y ha sido citado en más de 10 publicaciones. El anticuerpo CDX2 ha sido validado por Sigma utilizando RNAseq ortogonal.

Acumulación de edad e inferencia basada en IHC de REG4

Se utilizó el paquete RHClones (https://github.com/ElEd2/RHClones) para estimar la pendiente relacionada con la edad de los clones detectados (N = 291). El mismo paquete se utilizó para inferir la dinámica de los clones basándose en los datos de IHC de REG4.

LCM para muestras positivas para REG4

Se cortaron secciones seriadas de muestras con regiones positivas para REG4 identificadas en secciones de LCM, como se describió anteriormente (sección de Micro-seq). Las regiones positivas para REG4, junto con algunos criptas WT adyacentes para aumentar la cantidad de material genómico, se cortaron y se secuenciaron utilizando la secuenciación de próxima generación para las mutaciones de KRAS (G12/G13, A146, Q61) utilizando los cebadores mencionados anteriormente (Tabla suplementaria 3). Se incluyeron muestras WT/negativas para REG4 para cada paciente para excluir cualquier posible SNP. Para la llamada de mutaciones, se utilizó un ensamblador de pares finales para secuencias de Illumina, PANDAseq, para fusionar los pares de lecturas correspondientes en un amplicon in silico^[42]. Se identificaron las diez especies diferentes principales que comenzaron y terminaron con la secuencia de cebador correcta y se calculó su frecuencia. Las muestras mutadas contenían especies WT junto con especies KRAS mutadas.

Procesamiento de la diapositiva Visium y preparación de la biblioteca

Se cortaron secciones seriadas con los clones mutados de interés de MicroSeq y AceStripseq, con la ROI colocada en el centro de la diapositiva de vidrio. El procesamiento de la diapositiva y la preparación de la biblioteca se realizaron de acuerdo con el protocolo 10X Visium Cytassist FFPE.

Análisis de datos de Visium

La alineación de datos se realizó utilizando Spaceranger (10X). Se utilizó LoupeBrowser para dibujar alrededor de las regiones KRAS mutadas y las regiones WT coincidentes y extraer los códigos de barras. Los puntos con UM1 < 1000 genes < 900 y porcentaje de genes mitocondriales > 20% se eliminaron del análisis. El análisis de componentes principales, el agrupamiento no supervisado, el análisis de expresión génica diferencial y la reducción de dimensionalidad se realizaron utilizando la biblioteca de experimentos espaciales en R.

La firma de mutación de KRAS se generó compilando los 400 genes con mayor expresión diferencial de cada uno de los 14 clones y seleccionando los genes comunes a al menos dos clones (N = 628 genes únicos). Se realizó GSEA utilizando las bibliotecas Cluster Profiler y msigdbr en R.

Dado que identificamos que R465C se asocia con un fuerte sesgo de selección, utilizamos la información de localización espacial de los clones proporcionada por Micro-seq y AceStrip-seq para realizar un análisis de reducción de dimensionalidad en las regiones más pequeñas que contienen clones mutantes de FBXW7 R465C. No apareció ningún grupo distinto en dichas regiones, lo que probablemente refleja la falta de alteraciones transcripcionales asociadas con la pérdida de un alelo supresor tumoral y/o el pequeño tamaño de los clones, lo que impidió la agrupación.

Pruebas estadísticas

Todas las pruebas estadísticas se realizaron utilizando las funciones predeterminadas disponibles en el lenguaje de programación R. Se realizó la prueba exacta de Fisher utilizando la función fisher.exact() con los valores de parámetros predeterminados.

Resumen de la información presentada

Se puede encontrar información adicional sobre el diseño de la investigación en el Resumen de información de Nature Portfolio vinculado a este artículo.