Un conjunto de métodos de fusión de datos multimodales basados en aprendizaje automático para análisis de riesgo pronóstico en cáncer colorrectal metastásico tratado con bevacizumab.

El avance de la medicina de precisión se ha visto significativamente acelerado por los últimos avances en bioinformática y tecnologías genómicas, lo que exige la integración de datos multiómicos y clínicos (multimodales) utilizando métodos de aprendizaje automático de vanguardia para la investigación oncológica traslacional. Este enfoque transformador se centra en el aprovechamiento de técnicas de aprendizaje automático en el entorno clínico para descubrir biomarcadores multimodales, con el objetivo final de prolongar la vida de los pacientes y mejorar su calidad de vida. Recientemente, varios estudios innovadores han demostrado el potencial del aprendizaje automático y la inteligencia artificial (IA) en el cáncer colorrectal (CCR). Por ejemplo, Liu et al. (2025)^[1] han empleado un enfoque de descubrimiento de fármacos multi-objetivo basado en el aprendizaje automático, utilizando firmas de biomarcadores para identificar objetivos terapéuticos de precisión en el cáncer de colon, lo que ejemplifica cómo el aprendizaje automático puede revelar patrones moleculares relevantes. Además, Bahrambanan et al. (2025)^[2] desarrollaron un sistema de clasificación basado en el aprendizaje profundo para predecir la respuesta a la radiocuemioterapia en el CCR, demostrando que las arquitecturas avanzadas de IA superan al aprendizaje automático convencional en la modelización de la sensibilidad al tratamiento. De manera similar, Li et al. (2025)^[3] aplicaron la modelización cuantitativa de sistemas del ciclo de inmunidad contra el cáncer para predecir la progresión de la enfermedad en el CCR metastásico (m) avanzado, demostrando cómo los marcos basados en el aprendizaje automático y la modelización mecanicista pueden capturar la dinámica tumoral-inmune relevante para los resultados clínicos. Aquí, hemos investigado la aplicación de nuestra herramienta de aprendizaje automático, Phenotype Mapping (PhenMap)^[4], para identificar candidatos de biomarcadores multimodales en pacientes con CCR, seguido de la puntuación de riesgo pronóstico basada en el aprendizaje automático para asociar los biomarcadores identificados con la resistencia a la terapia combinada con bevacizumab (BVZ).

El CCR se clasifica como el tercer cáncer más diagnosticado en hombres y mujeres y está relacionado con una mortalidad y morbilidad significativas^[5],^[6]. Casi el 50% de los pacientes diagnosticados desarrollarán enfermedad metastásica. El enfoque de tratamiento actual para el mCCR mutante en RAS implica la quimioterapia estándar de atención (SOC) basada en 5-fluorouracilo combinada con el inhibidor de la angiogénesis BVZ. Los resultados de los ensayos clínicos de fase III han demostrado que la adición de BVZ a la quimioterapia mejora las tasas de respuesta y prolonga la supervivencia de los pacientes con mCCR^[7]–^[9]. Sin embargo, la respuesta solo se observa en un subconjunto de pacientes, y el beneficio clínico general de BVZ es limitado, lo que conduce a una eventual progresión de la enfermedad^[8]. La terapia con BVZ también se asocia con toxicidades no deseadas y altos costos de tratamiento. No obstante, en ausencia de otros tratamientos, se administra BVZ a los pacientes con mCCR mutante en RAS como tratamiento de primera línea^[9],^[10]. Es importante destacar que, si bien nosotros y otros hemos propuesto previamente varias entidades genómicas como posibles predictores de la respuesta a BVZ^[11]–^[17], no existen biomarcadores disponibles para identificar a los pacientes con mCCR que no responderán a BVZ. Para abordar este problema clínico pendiente, aquí hemos desarrollado una secuencia de trabajo de selección de características y metodologías de puntuación de riesgo basadas en el aprendizaje automático para detectar tanto a los pacientes que responden a BVZ como a los que no. Específicamente, hemos aprovechado nuestra herramienta de aprendizaje automático no supervisado previamente establecida, PhenMap^[4], que admite la integración simultánea de datos multiómicos y fenómicos del cáncer. Esta herramienta puede facilitar aún más la identificación de subtipos clínicamente relevantes y funcionales asociados con puntos finales específicos. Por ejemplo, en el contexto del cáncer de mama, PhenMap ha identificado previamente subtipos para la predicción de la respuesta a un inhibidor de la quinasa dependiente de ciclinas (CDK) 4/6^[4].

Aquí, PhenMap se aplicó a la cohorte clínica ANGIOPREDICT^[18], una colección retrospectiva agrupada de tumores (y metadatos clínicos asociados) de pacientes que se sometieron principalmente a terapia combinada con BVZ, en centros de Irlanda, los Países Bajos y Alemania^[19]. En el presente estudio, aprovechamos conjuntos de datos de alteraciones del número de copias (CNA) y perfiles de mutación generados previamente para estos pacientes^[16],^[17],^[19], que se integraron con datos clínicos utilizando la herramienta de aprendizaje automático PhenMap. La salida de PhenMap se utilizó además para desarrollar grupos de riesgo pronóstico que identificaron biomarcadores asociados con pacientes que eran resistentes al tratamiento combinado con BVZ. Se muestra un diagrama de flujo completo que resume la secuencia de trabajo en la Fig. 1**. En general, describimos el desarrollo de un marco de aprendizaje automático que fusiona datos multimodales con el pronóstico del paciente antiangiogénico y las respuestas al tratamiento para predecir el resultado en el mCCR. Nuestro nuevo enfoque se alinea con la nueva era de la medicina de precisión contra el cáncer, que se ha visto reforzada por los recientes avances en la aplicación de la IA y las metodologías de bioinformática a los conjuntos de datos de la biología tumoral.

Métodos

Conjuntos de datos

Los datos generados a partir de la cohorte retrospectiva ANGIOPREDICT^[19] (tumores localmente irresecables o mCCR de pacientes tratados con BVZ y quimioterapia) se utilizaron en el presente estudio. En resumen, las muestras se habían sometido previamente a secuenciación del genoma completo de baja cobertura, y las puntuaciones de GISTIC 2.0^[20] para los datos de CNA generados^[17], donde se demostró que 107 regiones estaban significativamente delecionadas/amplificadas. Además, el estado de mutación génica se había derivado previamente de la secuenciación del exoma completo (ADN tumoral y germinal emparejado disponible)^[17]. Finalmente, también se recopilaron los datos clínicos asociados con cada paciente, incluidos los parámetros clinicopatológicos como el sexo, la edad, el estadio, la ubicación del tumor, la inestabilidad cromosómica (CIN), el número de puntos de ruptura y el número de mutaciones (consulte el estudio anterior^[17]). Estos conjuntos de datos (es decir, número de copias, mutación y clínicos) se filtraron para retener específicamente a los pacientes con información clínica, de mutación y de número de copias coincidentes. En total, se utilizaron datos de 117 pacientes con características coincidentes en los tres tipos de datos (107 regiones de CNA, 5 datos de mutación génica y datos clinicopatológicos) como datos de entrada para PhenMap.

PhenMap

PhenMap^[4] es un modelo de análisis factorial que modela datos observados de alta dimensión como una función lineal de variables de mapeo (VM) de baja dimensión correspondientes. El modelo intenta explicar esta variación capturada por las VM utilizando información fenotípica o covariables adicionales del paciente. Considere las mediciones, , que pueden representar la expresión de p genes del tumor i, y caracterizadas con l fenotipos, (aquí, 1 representa el intercepto), incluidos los parámetros clínicos y/o genotípicos, donde T representa la transpuesta de una matriz. PhenMap puede modelar el tumor i (i = 1… n tumores) correlacionado mediciones, yi, utilizando VM de baja dimensión, , donde q es el número de VM, asumiendo las siguientes distribuciones gaussianas; donde W es una matriz de cargas p x q que relaciona cada característica con q VM y es la matriz de covarianza residual diagonal. Los coeficientes de regresión, β, cuantifican el efecto de los fenotipos (covariables) en cada VM y es la matriz de covarianza residual de las VM. Este modelo se estima utilizando un marco bayesiano para abordar el problema del tamaño de muestra pequeño común en la genómica. Consulte Nyamundanda et al^[4] para obtener más detalles sobre la estimación de este modelo.

Para modelar la correlación entre las características, introducimos la siguiente prioridad de determinación automática de relevancia (ARD) en W, que impone una escasez a nivel de elemento para la selección de características, donde G es una distribución gamma, j = 1…pt, y k = 1…q. La prioridad ARD en el hiperparámetro de precisión λjk controla la contribución, wjk, de la característica j en la k-ésima VM. El aumento de λjk reduce wjk, la carga correspondiente, acercándola a cero, lo que induce una escasez en W.

Para permitir la selección de fenotipos (covariables), se estableció una prioridad g en los coeficientes de regresión para regularizar β y también para eliminar la sensibilidad de la estimación de parámetros a los cambios en las escalas de los fenotipos. Se consideró la siguiente prioridad g, donde g controla la covarianza previa.

Aplicación de PhenMap

Los datos clínicos, de mutación y de número de copias (detallados anteriormente) se integraron utilizando PhenMap, una herramienta de análisis factorial bayesiano disperso previamente descrita que modela simultáneamente los datos ómicos y fenómicos para la subtipificación funcional^[4]. Las puntuaciones de GISTIC del número de copias se trataron como características de salida, mientras que las mutaciones y las variables clínicas del paciente se trataron como covariables coincidentes para la muestra. La salida de PhenMap se utilizó además para desarrollar grupos de riesgo pronóstico que identificaron biomarcadores asociados con pacientes que eran resistentes al tratamiento combinado con BVZ. Se muestra un diagrama de flujo completo que resume la secuencia de trabajo en la Fig. 1**.

Para identificar los marcadores genómicos más sólidos que predicen la supervivencia libre de progresión (SLP) en pacientes con BVZ, las regiones de CNA y las mutaciones asociadas con las VM pronósticas se seleccionaron y procesaron mediante la regresión de Cox penalizada (basada en la red elástica)^[21]. Esto se realizó dividiendo repetidamente los datos en datos de entrenamiento/prueba y ajustando la regresión de Cox penalizada 100 veces para eliminar el sesgo de muestreo. Las características (regiones de CNA y mutaciones) seleccionadas por la regresión de Cox penalizada (basada en la red elástica) con una frecuencia del 100% se eligieron como posibles marcadores pronósticos. Para identificar los atributos genómicos que son altamente predictivos del pronóstico, ajustamos un modelo de regresión de Cox multivariada estándar con los marcadores seleccionados y seleccionamos aquellos que fueron significativos a un nivel de significación del 5% (eliminando iterativamente las características no significativas hasta que todas las características fueran significativas). Finalmente, la puntuación de riesgo (riesgo de mortalidad), Ci, para el i-ésimo paciente que recibió BVZ se calculó de la siguiente manera; donde g es el número de marcadores genómicos seleccionados, es el coeficiente de regresión (logaritmo natural de la razón de riesgos [HR]) para los valores del marcador j. Por lo tanto, las puntuaciones de riesgo son un producto de los marcadores genómicos seleccionados con su correspondiente HR (que representa el riesgo asociado con cada marcador).

Para asignar a los pacientes a grupos de riesgo, las puntuaciones de riesgo se dividieron en cuantiles (0,1 y 0,9) para separar el 10% inferior y superior de las puntuaciones, respectivamente, en las categorías de bajo riesgo y alto riesgo. A menudo, los extremos o las colas de estas distribuciones contienen muestras atípicas que se comportan de manera diferente a las demás, como los pacientes de bajo y alto riesgo. Los pacientes con puntuaciones intermedias se asignaron al grupo de riesgo moderado. La asociación de los grupos de riesgo con el resultado clínico (SLP) se evaluó utilizando modelos de Cox multivariados de riesgos proporcionales, junto con la edad, la ubicación del tumor, el estadio, el CIN y la inclusión de nuestros subtipos ANGIOPREDICT previamente publicados (derivados solo de las CNA)^[17].

Para los análisis de supervivencia, se utilizó el método de Kaplan-Meier para trazar las curvas de supervivencia. El punto final estudiado fue la SLP. La SLP fue el tiempo desde el inicio del estudio hasta la fecha de la primera recurrencia o la fecha en que se realizó el último seguimiento. Todos los análisis estadísticos se realizaron en R Versión 4.3.0.

Conjuntos de datos

Los datos generados a partir del estudio retrospectivo ANGIOPREDICT^[19] (tumores localmente irresecables o mCRC de pacientes tratados con BVZ y quimioterapia) se utilizaron en el presente estudio. En resumen, las muestras se habían sometido previamente a una secuenciación del genoma completo de baja cobertura, y las puntuaciones de GISTIC 2.0^[20] para los datos de CNA generados^[17] se utilizaron para identificar las alteraciones del número de copias, donde se demostró que 107 regiones estaban significativamente delecionadas/amplificadas. Además, el estado mutacional de los genes se había derivado previamente de la secuenciación del exoma completo (ADN tumoral y de la línea germinal emparejados disponibles)^[17]. Finalmente, también se recopilaron los datos clínicos asociados con cada paciente, incluidos los datos clinicopatológicos como el sexo, la edad, el estadio, la ubicación del tumor, la inestabilidad cromosómica (CIN), el número de puntos de ruptura y el número de mutaciones (ver estudio previo^[17]). Estos conjuntos de datos (es decir, número de copias, mutación y datos clínicos) se filtraron para retener específicamente a los pacientes con información clínica, mutacional y de número de copias coincidentes. En total, se utilizaron los datos de 117 pacientes con características coincidentes en los tres tipos de datos (107 regiones de CNA, 5 datos de mutación de genes y datos clinicopatológicos) como datos de entrada para PhenMap.

PhenMap

PhenMap^[4] es un modelo de análisis factorial que modela los datos observados de alta dimensión como una función lineal de las variables de mapeo de baja dimensión correspondientes (MV). El modelo intenta explicar esta variación capturada por las MV utilizando información fenotípica o covariables adicionales del paciente. Considere las mediciones, , que pueden representar la expresión de p genes del tumor i, y caracterizadas con l fenotipos (aquí, 1 representa la intersección), incluidos los parámetros clínicos y/o genotípicos, donde T representa la transpuesta de una matriz. PhenMap puede modelar el tumor i (i = 1… n tumores) utilizando MV de baja dimensión, , donde q es el número de MV, asumiendo las siguientes distribuciones gaussianas; donde W es una matriz de cargas p x q que relaciona cada característica con q MV y es la matriz de covarianza residual diagonal. Los coeficientes de regresión, β, cuantifican el efecto de los fenotipos (covariables) en cada MV y es la covarianza residual de las MV. Este modelo se estima utilizando un marco bayesiano para abordar el problema del tamaño de muestra pequeño común en la genómica. Consulte Nyamundanda et al^[4] para obtener más detalles sobre la estimación de este modelo.

Para modelar la correlación entre las características, introducimos la siguiente determinación automática de relevancia (ARD) en la prioridad de dispersión en W, que impone la dispersión por elementos para la selección de características, donde G es una distribución gamma, j = 1…pt, y k = 1…q. La prioridad ARD en el hiperparámetro de precisión λjk controla la contribución, wjk, de la característica j en la k-ésima MV. El aumento de λjk reduce wjk, la carga correspondiente, acercándola a cero, lo que induce la dispersión en W.

Para permitir la selección de fenotipos (covariables), se estableció una prioridad g en los coeficientes de regresión para regularizar β y también para eliminar la sensibilidad de la estimación de parámetros a los cambios en las escalas de los fenotipos. Se consideró la siguiente prioridad g, , donde g controla la covarianza prioritaria.

Aplicación de PhenMap

Los datos clínicos, mutacionales y de número de copias (detallados anteriormente) se integraron utilizando PhenMap, una herramienta de aprendizaje automático de análisis factorial bayesiano disperso previamente descrita que modela simultáneamente los datos ómicos y fenómicos para la subtipificación funcional^[4]. Las puntuaciones GISTIC del número de copias se trataron como características de salida, mientras que las mutaciones y las variables clínicas del paciente se trataron como covariables coincidentes para la muestra. La información pronóstica asociada entre las MV de PhenMap y la supervivencia libre de progresión (SLP) para los pacientes se modeló utilizando el modelo estándar de Cox de riesgos proporcionales.

Selección de características significativas y puntuación de riesgo

En primer lugar, para identificar los marcadores genómicos más robustos que predicen la SLP en los pacientes tratados con BVZ, las regiones de CNA y las mutaciones asociadas con las MV pronósticas se seleccionaron y procesaron mediante la regresión de Cox penalizada (basada en la red elástica)^[21]. Esto se realizó dividiendo repetidamente los datos en datos de entrenamiento/prueba y ajustando la regresión de Cox penalizada 100 veces para eliminar el sesgo de muestreo. Las características (regiones de CNA y mutaciones) seleccionadas por la regresión de Cox penalizada (basada en la red elástica) con una frecuencia del 100% se eligieron como posibles marcadores pronósticos. Para identificar los atributos genómicos que son altamente predictivos del pronóstico, ajustamos un modelo de regresión de Cox multivariante estándar con los marcadores seleccionados y seleccionamos aquellos que fueron significativos a un nivel de significación del 5% (eliminando iterativamente las características no significativas hasta que todas las características fueran significativas). Finalmente, la puntuación de riesgo (riesgo de mortalidad), Ci, para el i-ésimo paciente que recibió BVZ se calculó de la siguiente manera; donde g es el número de marcadores genómicos seleccionados, es el coeficiente de regresión (logaritmo natural de la razón de riesgos [HR]) para los valores del marcador j. Por lo tanto, las puntuaciones de riesgo son un producto de los marcadores genómicos seleccionados con su correspondiente HR (que representa el riesgo asociado con cada marcador).

Para asignar a los pacientes a grupos de riesgo, las puntuaciones de riesgo se dividieron en cuantiles (0,1 y 0,9) para separar el 10% inferior y superior de las puntuaciones, respectivamente, en las categorías de bajo riesgo y alto riesgo. A menudo, los extremos o las colas de estas distribuciones contienen muestras atípicas que se comportan de manera diferente a las demás, como los pacientes de bajo y alto riesgo. Los pacientes con puntuaciones intermedias se asignaron al grupo de riesgo moderado. La asociación de los grupos de riesgo con el resultado clínico (SLP) se evaluó utilizando modelos de Cox multivariantes de riesgos proporcionales, junto con la edad, la ubicación del tumor, el estadio, la CIN e incluyendo nuestros subtipos ANGIOPREDICT previamente publicados (derivados solo de las CNA)^[17].

Para los análisis de supervivencia, se utilizó el método de Kaplan-Meier para trazar las curvas de supervivencia. El punto final estudiado fue la SLP. La SLP fue el tiempo desde el inicio del estudio hasta la fecha de la primera recurrencia o la fecha en que se realizó el último seguimiento. Todos los análisis estadísticos se realizaron en R Versión 4.3.0.

Resultados

PhenMap selecciona características multimodales integradas en muestras de mCRC tratadas con BVZ

La integración conjunta basada en PhenMap de las CNA (como resultados) y las mutaciones y los datos demográficos/resultados clínicos (como covariables) (datos originales derivados del estudio ANGIOPREDICT^[17]) reveló las diez MV superiores (MV1 a MV10) (consulte la sección de métodos para obtener una explicación de las MV) (Fig. 2a) que representaron distintos contextos de heterogeneidad molecular en la cohorte ANGIOPREDICT entre las CNA y los datos de mutación.

Curiosamente, MV5 y MV8 se asociaron con la SLP (pronóstico; según el análisis de Cox de riesgos proporcionales de la SLP) (Fig. 2a). Para derivar de manera eficiente interpretaciones clínicas importantes, PhenMap identificó características integradas, CNA y covariables (fenotipos), responsables de impulsar la heterogeneidad en MV5 y MV8 (Fig. 2b y c). Las principales regiones de CNA significativamente asociadas con MV5 fueron 15q26.1, 2q33.1, 15q22.33, 15q11.2 y 15q21.1 (a un nivel de significación del 0,25%). Además, se descubrió que BRAF fue la única covariable significativa que impulsó la heterogeneidad de MV5 (a un nivel de significación del 5%). De manera similar, las principales regiones de CNA significativamente asociadas con MV8 (a un nivel de significación del 0,25%) fueron 1p36.31, 4q13.3, 1p36.11, 4p16.2, 1p33, 1p31.1, 1p34.2 y 1p21.1. No se encontraron mutaciones o covariables clínicas asociadas con esta MV (Fig. 2b y c). Por lo tanto, PhenMap seleccionó características pronósticas significativas asociadas con MV5 y MV8 que se utilizaron para el siguiente conjunto de análisis posteriores.

La firma de CNA/mutaciones identifica grupos de riesgo en mCRC tratado con BVZ

Para derivar un efecto conjunto de estas relaciones (entre CNA y mutaciones) sobre el pronóstico, las CNA y las mutaciones integradas se redujeron a una puntuación de riesgo, es decir, ponderando el riesgo de mortalidad para cada paciente en función de la SLP después del tratamiento con la terapia combinada de BVZ. Esta reducción de características se realizó primero ajustando la regresión de Cox penalizada (con un enfoque de aprendizaje automático de red elástica^[22]) en todas las CNA/mutaciones seleccionadas por PhenMap asociadas con MV5 y MV8 en 117 pacientes que recibieron BVZ. Se encontró que las seis CNA/mutaciones estaban asociadas (frecuencia del 100%) con la SLP en el modelo de regresión de Cox de red elástica (Fig. 2d). El análisis de regresión de Cox multivariante estándar sobre las seis CNA/mutaciones seleccionadas redujo el número de características relevantes a dos regiones de CNA (deleciones de 15q21.1 y 1p36.31) y el estado de mutación de BRAF. Las HR de las tres características finales se muestran en la Fig. 2e.

Se generó una puntuación de riesgo para cada paciente multiplicando los valores normalizados de las puntuaciones GISTIC de las dos CNA y la mutación de BRAF por sus correspondientes HR. La relación entre las puntuaciones de riesgo y sus MV de origen (MV5 y MV8) se muestra como un diagrama de dispersión con coeficientes de correlación y valores p en la Fig. 3a y b, mientras que la Fig. 3c visualiza estas asociaciones en un formato combinado, lo que ilustra la derivación de la firma de estas dos MV pronósticas. Posteriormente, todos los pacientes se dividieron en tres grupos de riesgo (bajo, moderado y alto; grupos designados CRC-BVZ-RiskAssigner) en función de sus puntuaciones de riesgo mediante cuantiles, asignando el 10% inferior de las puntuaciones de riesgo al grupo de bajo riesgo y el 10% superior de las puntuaciones de riesgo al grupo de alto riesgo, con el resto de los pacientes constituyendo el grupo de riesgo moderado (Fig. 4a). El riesgo de mortalidad para el grupo de alto riesgo (n = 12) fue mayor que el del grupo de bajo riesgo (n = 12) (p = 1,29e-07 HR: 10,78) (Fig. 4b). Además, hubo una diferencia pronóstica significativa entre los grupos de bajo riesgo y de riesgo moderado (p = 0,034 HR: 1,99). En general, el rendimiento del modelo de puntuación CRC-BVZ-RiskAssigner (resultante del análisis de las características del receptor operativo (ROC)) mostró un área bajo la curva (AUC) de 0,685 (IC del 95%: 0,672-0,698) (Fig. 4c).

En un análisis multivariante que incorporó la edad, la ubicación del tumor, el estadio, la CIN e incluyendo los subtipos ANGIOPREDICT^[17], los grupos de riesgo demostraron ser un predictor independiente de la SLP en los pacientes tratados con BVZ (Fig. 5a). El índice de concordancia (C-índice) del grupo CRC-BVZ-RiskAssigner es de 0,60 (IC del 95%: 0,56, 0,65), mientras que los grupos basados solo en los subtipos ANGIOPREDICT tienen un C-índice de 0,58 (Fig. 5b). Los intervalos de confianza de la CIN, el estadio y la ubicación no fueron mayores de 0,5 del C-índice (no significativo). En esencia, hemos demostrado que el grupo de alto riesgo solo contiene pacientes que no responden a BVZ cuando la respuesta del paciente se asoció con los grupos de riesgo. El grupo de bajo riesgo contenía a la mayoría de los respondedores (10 de 12 pacientes; 83,33%). El grupo de riesgo moderado comprendía tanto respondedores como no respondedores (Fig. 5c).

Curiosamente, MV5 y MV8 se asociaron con la supervivencia libre de progresión (SLP) (pronóstico; basado en el análisis de Cox de la SLP) (Fig. 2a). Para derivar de manera eficiente interpretaciones clínicas importantes, PhenMap identificó características integradas, alteraciones del número de copias (ANC) y covariables (fenotipos) responsables de impulsar la heterogeneidad en MV5 y MV8 (Fig. 2b y c). Las principales regiones de ANC significativamente asociadas con MV5 fueron 15q26.1, 2q33.1, 15q22.33, 15q11.2 y 15q21.1 (a un nivel de significancia del 0.25%). Además, se encontró que BRAF fue la única covariable significativa que impulsó la heterogeneidad de MV5 (a un nivel de significancia del 5%). De manera similar, las principales regiones de ANC significativamente asociadas con MV8 (a un nivel de significancia del 0.25%) fueron 1p36.31, 4q13.3, 1p36.11, 4p16.2, 1p33, 1p31.1, 1p34.2 y 1p21.1. No se encontraron mutaciones o covariables clínicas asociadas con este MV (Fig. 2b y c). Por lo tanto, PhenMap seleccionó características pronósticas significativas asociadas con MV5 y MV8 que se utilizaron para el siguiente conjunto de análisis posteriores.

La firma de ANC/mutaciones identifica grupos de riesgo en pacientes con mCRC tratados con BVZ

Para derivar un efecto conjunto de estas relaciones (entre ANC y mutaciones) sobre el pronóstico, las ANC y las mutaciones integradas se redujeron a una puntuación de riesgo, es decir, ponderando el riesgo de mortalidad para cada paciente en función de la SLP después del tratamiento con terapia combinada con BVZ. Esta reducción de características se realizó primero ajustando una regresión de Cox penalizada (con un enfoque de aprendizaje automático de red elástica ^[22]) en todas las ANC/mutaciones seleccionadas por PhenMap asociadas con MV5 y MV8 en 117 pacientes que recibieron BVZ. Seis ANC/mutaciones estuvieron asociadas (frecuencia del 100%) con la SLP en el modelo de regresión de red elástica (Fig. 2d). El análisis de regresión de Cox multivariante estándar sobre las seis ANC/mutaciones seleccionadas redujo el número de características relevantes a dos regiones de ANC (deleciones en 15q21.1 y 1p36.31) y el estado mutacional de BRAF. Los riesgos relativos (RR) de las tres características finales se muestran en la Fig. 2e.

Se generó una puntuación de riesgo para cada paciente multiplicando los valores normalizados de las puntuaciones GISTIC para las dos ANC y la mutación de BRAF por sus respectivos RR. La relación entre las puntuaciones de riesgo y sus MVs de origen (MV5 y MV8) se muestra como un diagrama de dispersión con coeficientes de correlación y valores p en la Fig. 3a y b, mientras que la Fig. 3c visualiza estas asociaciones en un formato combinado, lo que ilustra la derivación de la firma a partir de estos dos MVs pronósticos. Posteriormente, todos los pacientes se dividieron en tres grupos de riesgo (bajo, moderado y alto; designados como grupos CRC-BVZ-RiskAssigner) en función de sus puntuaciones de riesgo mediante cuantiles, asignando el 10% más bajo de las puntuaciones de riesgo al grupo de bajo riesgo y el 10% más alto de las puntuaciones de riesgo al grupo de alto riesgo, con el resto de los pacientes constituyendo el grupo de riesgo moderado (Fig. 4a). El riesgo de mortalidad para el grupo de alto riesgo (n = 12) fue mayor que el del grupo de bajo riesgo (n = 12) (p = 1.29e-07; RR: 10.78) (Fig. 4b). Además, hubo una diferencia pronóstica significativa entre los grupos de bajo y moderado riesgo (p = 0.034; RR: 1.99). En general, el rendimiento del modelo de puntuación CRC-BVZ-RiskAssigner (resultante del análisis de curvas ROC) mostró un área bajo la curva (AUC) de 0.685 (IC del 95%: 0.672–0.698) (Fig. 4c).

En un análisis multivariante que incorporó la edad, la ubicación del tumor, el estadio, la inestabilidad genómica (IG) e incluyó los subtipos ANGIOPREDICT ^[17], los grupos de riesgo demostraron ser un predictor independiente de la SLP en pacientes tratados con BVZ (Fig. 5a). El índice de concordancia (C-índice) del grupo CRC-BVZ-RiskAssigner es de 0.60 (IC del 95%: 0.56, 0.65), mientras que los grupos basados únicamente en los subtipos ANGIOPREDICT tienen un C-índice de 0.58 (Fig. 5b). Los intervalos de confianza de la IG, el estadio y la ubicación no fueron mayores de 0.5 (no significativo). Fundamentalmente, hemos demostrado que el grupo de alto riesgo solo contiene pacientes que no responden a BVZ cuando la respuesta del paciente se asoció con los grupos de riesgo. El grupo de bajo riesgo contenía a la mayoría de los respondedores (10 de 12 pacientes; 83.33%). El grupo de riesgo moderado estaba compuesto por respondedores y no respondedores (Fig. 5c).

Discusión

Aquí presentamos una completa plataforma de aprendizaje automático de múltiples pasos que fusiona conjuntos de datos de ómicas multimodales y resultados clínicos, junto con información pronóstica y respuesta al tratamiento para predecir el resultado del paciente después de la terapia combinada antiangiogénica en el entorno de mCRC. Específicamente, al combinar un método publicado previamente (PhenMap ^[4]) de integración de datos multimodales y fenómicos, con enfoques de aprendizaje automático de predicción de riesgos, identificamos grupos de riesgo pronósticos a partir de una cohorte de pacientes con mCRC que recibieron terapia combinada con BVZ [17, 19]. El riesgo de mortalidad fue significativamente mayor para los pacientes en el grupo de alto riesgo, con el 100% (n = 12) de los pacientes que no mostraron respuesta a BVZ, en comparación con el grupo de bajo riesgo (16.67%, 2 de 12). En general, la asignación al grupo de alto riesgo fue un predictor negativo independiente de la supervivencia en pacientes con mCRC tratados con BVZ.

La integración simultánea de datos tumorales multimodales y fenómicos, a través de la aplicación de PhenMap, identificó diez MVs, que asociaron características de diferentes tipos de datos de ómicas y covariables clínicas/mutacionales en diferentes contextos. Dos de estos MVs, a saber, MV5 y MV8, demostraron estar significativamente asociados con la SLP y, por lo tanto, se seleccionaron para el análisis posterior. Curiosamente, observamos que MV5 está negativamente asociado con el estado de BRAF, que se ha asociado con un mal pronóstico en mCRC ^[23]. De hecho, todos los pacientes en el grupo de alto riesgo presentaban la mutación tumoral de BRAF (junto con las otras dos ANC) y no respondieron a la terapia con BVZ. En otro lugar, datos previos del estudio clínico TRIBE (NCT00719797) indicaron que los pacientes con tumores mutados en BRAF tenían una mediana de supervivencia general (SG) más larga cuando fueron tratados con FOLFOXIRI + BVZ (N = 16 pacientes; SG de 19.0 meses) en comparación con FOLFIRI + BVZ (N = 12 pacientes; SG de 10.7 meses). Sin embargo, en general, estos pacientes tuvieron una SG más corta (con cualquiera de los regímenes) en comparación con los pacientes con RAS y BRAF de tipo salvaje (N = 45; SG de 33.5 meses [FOLFIRI + BVZ]), N = 48; SG de 41.7 meses [FOLFOXIRI + BVZ]) ^[24]. En este contexto, y particularmente cuando se combina con alteraciones genómicas adversas como las deleciones en 15q21.1 o 1p36.31, la mutación de BRAF sigue estando fuertemente asociada con un mal resultado y parece estar enriquecida en nuestro grupo de no respondedores de alto riesgo.

Para cada MV, PhenMap proporcionó asociaciones de variables fenotípicas con patrones de ANC. Aquí, se demostró que las alteraciones en el cromosoma 15q tenían un impacto significativo en la asociación con MV5, mientras que las alteraciones en el cromosoma 1p se asociaron con MV8. Curiosamente, se ha demostrado previamente que las deleciones en los brazos largos del cromosoma 15 (15q) son indicadores de mal pronóstico y progresión de la enfermedad en pacientes con mCRC [25, 26], lo que sugiere la presencia de importantes genes supresores de tumores. Específicamente, se sabe que la deleción en 15q21.1 está asociada con la pérdida del gen B2M, que codifica para el complejo mayor de histocompatibilidad (CMH-I) ^[27]. Además, Sheffer et al. proponen un papel de la deleción en 1q en las alteraciones de las vías de la fosforilación oxidativa que conducen a un peor pronóstico ^[25]. Se justifica una mayor investigación sobre el papel de la pérdida de 15q y 1q.

Después de reducir el número de características específicas de MV mediante la selección de aquellas que están más fuertemente asociadas con la SLP, utilizamos características significativas de MV5 y MV8 para calcular una "puntuación de riesgo". A continuación, diseñamos un sistema de clasificación de la puntuación de riesgo, en el que los pacientes se asignan a tres grupos de riesgo diferentes (bajo, moderado y alto) en función de su respuesta prevista a la terapia combinada con BVZ. Observamos que los pacientes con mCRC de alto y moderado riesgo demostraron una peor SLP en comparación con los pacientes con mCRC de bajo riesgo. Además, el grupo de alto riesgo demostró ser un predictor negativo independiente de la supervivencia en pacientes con mCRC tratados con BVZ. Por lo tanto, se pueden extraer las características asociadas con los pacientes de alto riesgo y potencialmente probarse como biomarcadores clínicos para el subtipo de alto riesgo. Esto permitiría a los clínicos identificar a los no respondedores y evitarles el tratamiento con BVZ y los efectos secundarios asociados, que afectan negativamente a la calidad de vida. En general, para desarrollar aún más nuestra novedosa plataforma metodológica, el biomarcador candidato y los grupos de riesgo ahora requieren una validación rigurosa en una cohorte separada y bien anotada de muestras de pacientes con mCRC. En última instancia, la herramienta PhenMap podría combinarse en el futuro con un ensayo de secuenciación de nueva generación (NGS), como el ensayo de diagnóstico complementario FoundationOne, que detecta tanto las ANC como la carga mutacional en una sola muestra del paciente. La salida de esta prueba aprobada por la FDA podría utilizarse como entrada para la plataforma metodológica PhenMap. Por lo tanto, al integrar datos de número de copias, mutaciones y respuesta clínica, se podrían acelerar las predicciones robustas del resultado del paciente.

Si bien los hallazgos de nuestro marco de aprendizaje automático integrado son prometedores, una serie de limitaciones merecen consideración. Primero, la cohorte ANGIOPREDICT, aunque está caracterizada de forma única, sigue siendo relativamente pequeña en tamaño de muestra para el modelado multi-ómico de alta dimensión; esto puede limitar la potencia estadística y podría hacer que ciertas asociaciones sean específicas de la cohorte. En segundo lugar, aunque PhenMap integra eficazmente los datos multimodales, el presente estudio se centró principalmente en las ANC, un conjunto limitado de mutaciones genéticas y las variables clínicas y patológicas básicas. Para capturar una biología adicional relevante para la resistencia a BVZ, puede ser necesario agregar capas como la transcriptómica, la epigenómica, la proteómica y las características del microentorno tumoral (por ejemplo, la patología espacial, las métricas de infiltración inmune). En tercer lugar, nuestro modelo predictivo se entrenó y se evaluó dentro de un único conjunto de datos retrospectivo y, por lo tanto, requiere validación en cohortes independientes recopiladas de forma prospectiva, idealmente incluyendo pacientes tratados bajo regímenes terapéuticos estandarizados y con flujos de trabajo de secuenciación armonizados. Por último, la interpretación mecanicista de la firma de ANC/mutaciones identificada, en particular las vías biológicas que vinculan la pérdida de 15q21.1 y 1p36.31 con la resistencia a la terapia antiangiogénica, sigue siendo incompleta. Los estudios futuros deberían priorizar la validación externa multicohorte, la disección experimental de las consecuencias mecanicistas de estas alteraciones y la exploración de cómo estos biomarcadores pueden interactuar con las combinaciones emergentes de terapias dirigidas e inmunoterapias. La expansión de la plataforma para integrar modalidades adicionales, incluyendo la patología digital y las características de la biopsia líquida, puede mejorar aún más la precisión predictiva y respaldar la traducción en herramientas de apoyo a la decisión clínicamente implementables.

En conclusión, hemos establecido un nuevo método de fusión basado en el aprendizaje automático para integrar datos de ómicas multimodales con resultados clínicos. Al emplear este método, hemos identificado un posible biomarcador combinado de ANC/mutaciones con puntuaciones de riesgo asociadas para la detección de pacientes con mCRC que no responderán a la terapia combinada con BVZ. Este posible biomarcador ahora requiere una mayor validación. Es importante destacar que la herramienta de aprendizaje automático presentada y la innovadora plataforma de análisis de fusión de datos multimodales ahora se pueden aplicar a otros conjuntos de datos multimodales en el entorno de la oncología.