Marco de trabajo multimodal para el análisis conjunto de datos de secuenciación de ARN de células individuales y de secuenciación del receptor de células T que predice la respuesta de las células T a la inmunoterapia contra el cáncer.

Las células T desempeñan un papel fundamental en la respuesta inmune adaptativa, y comprender y modular sus actividades se ha convertido en un objetivo central en la inmunoterapia^[1],^[2]. La reactividad de una célula T está determinada por su receptor de células T (TCR) único, que puede reconocer específicamente antígenos extraños presentados por moléculas MHC en la superficie celular de las células presentadoras de antígenos^[3],^[4]. En presencia de cáncer, las células T naïvas circulantes se activan en los ganglios linfáticos que drenan el tumor, lo que permite su migración hacia los tumores^[5]. Una célula T naïva se activa cuando reconoce un antígeno tumoral a través de su TCR y luego se somete a una expansión clonal y diferenciación, ya sea en una célula T efectora que mata las células cancerosas o en una célula T de memoria que protege durante exposiciones adicionales al mismo cáncer^[6]. Las células cancerosas adquieren muchos mecanismos para evadir la respuesta de las células T^[7]. Esta observación ha impulsado importantes avances en la inmunoterapia, incluido el bloqueo de puntos de control inmunitario, la terapia adoptiva de células T y las vacunas contra el cáncer, todos diseñados para restaurar la respuesta de las células T^[8]. Para desarrollar inmunoterapias contra el cáncer más eficaces, necesitamos una comprensión multifacética de la respuesta de las células T que tenga en cuenta que esta respuesta es a la vez muy dinámica y coordinada en múltiples sitios tisulares.

El advenimiento de la secuenciación de ARN de una sola célula (scRNAseq) ha mejorado enormemente la resolución de los subconjuntos de células T, revelando la heterogeneidad subyacente en los estados transcripcionales, especialmente en los tumores, con variación entre pacientes y entre los diferentes tipos de tumores. A pesar de esta heterogeneidad, los atlas transcriptómicos pancancerosos recientes que perfilan las células T infiltrantes tumorales en diferentes tipos de tumores han identificado estados y vías de desarrollo compartidos de las células T. Las señales del microambiente tumoral a menudo impulsan a las células T hacia un estado de agotamiento caracterizado por una mayor expresión de receptores inhibidores de células T (por ejemplo, PD1, CTLA4, LAG3), pérdida de la función de eliminación efectora, reducción de la secreción de citocinas y, en última instancia, evasión inmune tumoral^[9]–^[13]. Un tema recurrente es que los pacientes que no responden con éxito a la inmunoterapia contra el cáncer a menudo tienen células T infiltrantes tumorales en un estado de agotamiento^[14], lo que indica que los estados de las células T son un biomarcador predictivo importante.

Si bien la transcriptómica es potente para identificar los estados de las células T, la secuencia del TCR es, en última instancia, lo que impulsa el reconocimiento de los antígenos del cáncer y una respuesta inmune protectora. La secuenciación de TCR de una sola célula (scTCRseq) determina tanto las secuencias de la cadena α como de la cadena β del TCR que constituyen el TCR y nos permite definir la clonalidad en una población de células T, identificando los clones de células T con el mismo TCR que se han expandido en respuesta a un antígeno de clones que no se han expandido^[15],^[16]. Cuando se combina con scRNAseq de la misma célula, esto proporciona un vínculo potente entre un clonotipo de célula T y su estado celular, lo que ha dado lugar a muchas ideas sobre las actividades de los clones de células T expandidos. Por ejemplo, en un estudio de células T infiltrantes de glioma, el análisis de las células T expandidas clonalmente identificó un nuevo subconjunto que expresa el receptor inhibidor KLRB1, que, cuando se bloqueó, resultó en una mayor eliminación de las células de glioma^[17]. La secuencia del TCR no solo informa sobre el clonotipo de la célula T, sino que también puede actuar como un código de barras molecular para rastrear las células T a lo largo del tiempo y en diferentes sitios tisulares^[18]. Por ejemplo, en un ensayo clínico de bloqueo de puntos de control inmunitario neoadyuvante en carcinomas de células escamosas de cabeza y cuello (HNSCC), se recolectaron simultáneamente células T infiltrantes tumorales y células T de sangre periférica antes y después del tratamiento^[19]. Los análisis multimodales conjuntos de las secuencias de TCR y los transcriptomas de las células T identificaron células T CD8+ que se expandieron clonalmente durante la terapia, que habían existido previamente en el tumor antes del tratamiento, que expresaban programas de memoria y citotoxicidad residentes en los tejidos, y cuyos clonotipos se podían detectar sistémicamente en la sangre.

El perfilado multimodal del sistema inmunitario durante el tratamiento de la enfermedad y en diferentes sitios tisulares permite el descubrimiento de biomarcadores clínicos predictivos del pronóstico del paciente y la respuesta a la inmunoterapia. Un enfoque prometedor es vincular las mediciones en la sangre periférica con las mediciones en el tejido tumoral, ya que la recolección de sangre es un procedimiento clínico de rutina y se puede realizar repetidamente durante el curso del tratamiento del paciente^[20]. Por ejemplo, Epi proColon, aprobado por la FDA, analiza muestras de sangre en busca de ADN metilado de Septin 9 que se desprende de los tumores como biomarcador para el cáncer colorrectal^[21]. De manera similar, las células inmunitarias periféricas también pueden servir como biomarcadores del cáncer. Hasta la fecha, se han identificado múltiples biomarcadores de células T periféricas, incluido un vínculo entre la expansión clonal periférica e intratumoral en respuesta a la inmunoterapia anti-PD-1/PD-L1, la proliferación temprana de células T periféricas que predice la respuesta al bloqueo de puntos de control inmunitario y los marcadores de superficie celular NKG2D, CD39 y CX3CR1 en las células T periféricas, que identifican las células T infiltrantes tumorales en el melanoma^[22]–^[24].

Estos descubrimientos recientes en la capacidad predictiva de la sangre periférica con respecto a la dinámica tumoral subrayan la gran necesidad de integración y traducción de datos multimodales. Se perfilan rutinariamente diferentes modalidades de datos, incluidas las tabulares, el texto y las imágenes, en células y tejidos, y su integración proporciona información al unificar múltiples perspectivas biológicas^[25],^[26]. Sin embargo, se han realizado muy pocos estudios multimodales que recopilen muestras inmunitarias longitudinalmente durante el tratamiento del cáncer y en múltiples sitios tisulares, ya que son costosos, difíciles en el entorno clínico y requieren análisis integrativos novedosos^[19],^[27],^[28]. La integración de conjuntos de datos de diferentes modalidades se complica por la variabilidad técnica y las características únicas de cada modalidad. En el contexto de una sola célula, la heterogeneidad de múltiples distribuciones complejas presenta desafíos, pero también oportunidades para este tipo de análisis^[29]–^[31]. Por ejemplo, varios métodos, incluidos CoNGA, tessa y mvTCR, integran conjuntamente transcriptomas y secuencias de TCR de células individuales^[32]–^[34]. Además, un número creciente de métodos permite la traducción de una modalidad a otra^[35], por ejemplo, la predicción de datos de scRNAseq en células T CD4+ a partir de la accesibilidad de la cromatina de las células emparejadas^[29]. Muchos de estos métodos se basan en autoencoders variacionales (VAE), arquitecturas de redes neuronales eficaces para la integración multimodal que ofrecen una gran flexibilidad en el modelado de diferentes distribuciones de datos y su integración en una representación compartida^[36],^[37]. Sin embargo, hasta donde sabemos, no existen modelos multimodales para las células T que puedan predecir la respuesta tumoral a partir de las mediciones de la sangre periférica.

Presentamos un método de integración de datos multimodal que modela la scRNAseq y la scTCRseq emparejadas medidas en múltiples puntos temporales y en múltiples sitios tisulares (Fig. 1a). Al integrar estas dos modalidades de datos a través de una arquitectura VAE, podemos dilucidar la relación entre el transcriptoma y la clonalidad de las células T e identificar subconjuntos de células T funcionalmente relacionados (Fig. 1b). Nuestro modelo de integración de TCR y ARN (TRIM) no solo integra estas dos modalidades, sino que también puede generar predicciones para un punto temporal o sitio tisular no medido. Probamos TRIM en una cohorte clínica ampliamente perfilada de pacientes con HNSCC que reciben terapia de bloqueo de puntos de control inmunitario, con mediciones emparejadas de scRNAseq y scTCRseq del tumor y la sangre periférica antes y después del tratamiento^[19]. También probamos TRIM en un conjunto de datos similar que contiene datos emparejados de secuenciación de ARN y TCR de una sola célula de pacientes con cáncer colorrectal (CRC)^[27], y en un conjunto de datos de pan-cáncer emparejado de ARN y TCR recopilado de tumores y tejidos adyacentes^[11], lo que demuestra que nuestro modelo se generaliza más allá de los conjuntos de datos en los que se desarrolló. Nuestro marco predice con precisión cómo cambian las células T infiltrantes tumorales después del tratamiento analizando solo las células T de la sangre periférica antes del tratamiento. Esperamos que TRIM sea un marco de modelado potente en cohortes clínicas, donde el perfilado de todos los tejidos de interés en cada punto temporal es un desafío, así como en el descubrimiento de biomarcadores para el pronóstico del tumor y la respuesta al tratamiento.

Resultados

TRIM utiliza un autoencoder variacional condicional para predecir el efecto del tratamiento en el tumor a partir de la sangre pretratamiento

Para integrar la información de las secuencias de ARN y TCR en las células T y hacer predicciones sobre los estados de las células T en sitios y puntos temporales no vistos, proponemos TRIM, un modelo generativo basado en un VAE condicional (Fig. 1b, Métodos: TRIM). TRIM se entrena con datos de secuenciación de ARN y TCR de una sola célula emparejados de células T en la sangre circulante antes del tratamiento, y luego predice el estado de la población de células T a nivel de una sola célula tanto en los tumores pretratamiento como en los postratamiento, así como en las muestras de sangre postratamiento (Fig. 1b). TRIM aprende una representación conjunta de la expresión génica y la clonalidad utilizando dos codificadores separados para mapear el perfil de ARN y la secuencia de TCR de una célula en un espacio latente conjunto, condicionado a la fuente (sangre o tumor), el punto temporal (pretratamiento o postratamiento) y el paciente. Para garantizar que el espacio latente conjunto capture información de ambos, el perfil de ARN de la célula y su secuencia de TCR, TRIM utiliza dos decodificadores separados para reconstruir las dos modalidades (Fig. 2a).

Cuando están presentes múltiples modalidades emparejadas, los autoencoders tradicionales minimizan el error al traducir entre modalidades utilizando el codificador de una y el decodificador de la otra. Esta estrategia ha sido eficaz para integrar varias modalidades de una sola célula, como la secuenciación de ARN, la secuenciación ATAC y la imagen de la cromatina^[35],^[38],^[39]. Se aplica una pérdida de reconstrucción para alentar a que la salida se acerque a los puntos originales. Para el espacio de ARN, esta pérdida se calcula utilizando la norma L2 (Métodos). Las representaciones de TCR utilizadas en nuestro modelo se derivan de incrustaciones aprendidas a través de un autoencoder convolucional, que entrenamos como un paso de preprocesamiento. Específicamente, cada TCR se representa mediante la secuencia de aminoácidos CDR3 de su cadena β, codificada como una secuencia con relleno nulo donde cada aminoácido se representa mediante sus factores de Atchley, descriptores numéricos informados bioquímicamente que capturan propiedades como la polaridad, la carga y la hidrofobicidad^[40] (Métodos). Extraemos las incrustaciones de la capa latente media del autoencoder convolucional entrenado para que sirvan como entrada para el modelo TRIM. A pesar de centrarse en la secuencia CDR3 altamente variable de las cadenas β de los TCR, el modelo aprende implícitamente los segmentos génicos V y J contribuyentes, lo que da como resultado un espacio de representación aprendido organizado de manera biológicamente significativa (Fig. S1). Nos centramos en la cadena β por razones prácticas y biológicas. Primero, en la secuenciación V(D)J de una sola célula, el emparejamiento α- β a menudo es incompleto. En las muestras de sangre perfiladas en ref. ^[19], el 87% de las células tienen cadenas emparejadas, con β recuperada con más frecuencia que α (95% frente a 89%). En términos más generales, los grandes recursos públicos de TCR son predominantemente solo β, con el 94% de las secuencias que contienen solo la cadena β^[41]. En segundo lugar, los estudios anteriores que incluyen ambas cadenas informan solo ganancias de rendimiento marginales, probablemente porque la secuencia CDR3β es el determinante dominante de la especificidad péptido-MHC^[42],^[43]. Dadas estas limitaciones experimentales y la disponibilidad de datos, utilizamos la cadena β para maximizar la cobertura de la muestra y, al mismo tiempo, conservar la mayor parte de la señal relevante.

También observamos que los perfiles de secuenciación de ARN (RNA-seq) contienen información limitada sobre la secuencia del TCR de una célula, más allá de su clonalidad (como discutiremos con más detalle más adelante); por lo tanto, no es posible establecer una traducción directa entre el TCR y el ARN. Para abordar esta limitación, modelamos la clonalidad de las células T en lugar de las secuencias de TCR. Para preservar la estructura local en el espacio del TCR sin penalizar las variaciones no relacionadas, introducimos la pérdida de clonalidad contrastiva. Esta pérdida fomenta que las células T del mismo clon estén cerca unas de otras en el espacio latente, al tiempo que garantiza que diferentes clones de células T permanezcan separados por un margen (Métodos, pérdida de clonalidad contrastiva). Este enfoque desplaza el foco de la representación exacta del TCR a la preservación de la clonalidad.

Además, nuestro modelo incorpora las condiciones del paciente, la fuente y el punto temporal como entradas, representadas por vectores de longitud fija que se aprenden conjuntamente con el resto del modelo a través de una red de incrustación separada, F (Fig. 2a). En contraste con el enfoque estándar, que aprende vectores condicionales indexando en una matriz de incrustación inicializada aleatoriamente (con los parámetros de esta matriz aprendidos durante el entrenamiento), nuestro método está diseñado para aprender las representaciones vectoriales de las condiciones utilizando los perfiles de ARN y TCR asociados con cada condición. En concreto, la red F procesa un subconjunto de células de entrenamiento, Semb, y las mapea en un espacio latente de demb dimensiones. Para calcular el vector de incrustación para una condición específica, TRIM promedia las representaciones de las células correspondientes a esa condición. Por ejemplo, el vector del punto temporal pretratamiento se calcula promediando las representaciones de las células pretratamiento en Semb (Métodos, Aprendizaje de incrustaciones de condiciones). Para las incrustaciones de pacientes, nos centramos específicamente en las células sanguíneas pretratamiento, ya que estos son los datos disponibles para los pacientes no vistos. Las incrustaciones se aprenden a partir de los datos de ARN y TCR, utilizando dos redes, FRNA y FTCR, cuyas salidas se promedian para formar la incrustación final. Es importante destacar que nuestro modelo puede aprovechar el hecho de que la distribución de ARN o TCR de un paciente de una fuente/punto temporal es muy informativa sobre sus distribuciones en otros entornos. Demostramos que este enfoque es más eficaz para generar distribuciones precisas de ARN y TCR en diferentes condiciones y, lo que es crucial, permite una extrapolación precisa a nuevos pacientes no vistos. Durante la evaluación del modelo, discretizamos las salidas continuas del decodificador para estimar los clones de TCR, asignando pseudoclones utilizando un algoritmo voraz que agrupa los puntos dentro de un margen en el mismo clon (Fig. 2b, Métodos, Asignación de pseudoclones). Para maximizar la utilidad de los datos y probar su capacidad de generalización, entrenamos con todos los datos excepto los de un paciente y probamos el modelo con el paciente excluido. Repetimos este procedimiento en un proceso de validación cruzada para cada paciente, con el fin de obtener el conjunto completo de predicciones. Esto permite comparar la salida con la verdad fundamental en todo el conjunto de datos, pero todo ello proviene de predicciones fuera de la muestra.

Las secuencias de ARN y TCR de una sola célula están relacionadas solo a través del recuento de clones

Para explorar la relación entre el TCR y el ARN, analizamos datos de secuenciación de ARN de una sola célula (scRNAseq) y secuenciación de TCR de una sola célula (scTCRseq) emparejados de 27 pacientes con CHCC^[19]. El conjunto de datos incluye mediciones tanto de sangre periférica como de células T infiltrantes tumorales tomadas antes y después de la terapia con inhibidores de puntos de control inmunitario. Consta de 73.961 células T con expresión de ARN medida en 29.750 genes, emparejadas con secuencias de la cadena β de CDR3 del TCR. Si bien cada perfil de expresión de ARN es único, muchas células son clones con secuencias de TCR idénticas, y solo el 69,8 % de las células T tienen TCR únicos. Si bien estudios previos sugieren que los clones de células T idénticos comparten expresiones de ARN similares^[32], quedan muchas preguntas sobre la relación entre el espacio de expresión de ARN y el espacio de TCR.

Para ello, primero visualizamos los datos de scRNAseq utilizando UMAP^[44] (Fig. 3a). Se puede observar un fuerte solapamiento en las distribuciones de los pacientes individuales (Fig. 3a^[44], esquina inferior izquierda), con cuatro grupos distintos: las células T de la sangre están separadas de las del tumor, y las células T CD4+ son distintas de las células T CD8+ en ambas fuentes. Si bien los tipos de células están claramente separados en el espacio transcriptómico, el efecto del tratamiento es más sutil. No se observa un cambio discernible en las células sanguíneas, pero se produce un ligero cambio en las células tumorales, especialmente en las células CD8+. Al asignar a cada célula su recuento de clones (que corresponde al número de células que expresan el mismo TCR), se puede observar una mayor clonalidad en las regiones del espacio de expresión de ARN que corresponden a partes de los grupos de células T CD8+. Los tamaños de los clones más grandes son más frecuentes en las células CD8+ y en el tumor en comparación con las células CD4+ y las de la sangre (Fig. S2).

Para evaluar la hipótesis de que los clones están más cerca en la expresión de ARN, comparamos las células por distancia por pares en el espacio de ARN frente a la distancia de Hamming por pares en el espacio de TCR (Fig. 3b). La distribución de las distancias de ARN entre los clones de TCR idénticos (distancia de Hamming de cero) es sustancialmente menor que las distancias de ARN para las células que están más alejadas en el espacio de TCR. Curiosamente, una vez que superamos los clones de TCR idénticos, existe muy poca relación entre la secuencia de TCR y la expresión de ARN de dos células. Para encontrar alguna relación, debemos filtrar los clones con un recuento mayor que uno, e incluso entonces, la diferencia entre la distancia media de ARN entre los clones que están a 0 de Hamming y los que están a 1 de Hamming es 4 veces mayor que la diferencia entre las medias de las células que están a 1 de Hamming y a 10 de Hamming. En otras palabras, dos células cuyas secuencias de TCR difieren en una sola posición están tan lejos en el espacio de ARN como dos células con secuencias de TCR completamente diferentes. Observamos que esto es cierto no solo a nivel del conjunto de datos, sino también si se restringe solo a un paciente, a una fuente (sangre/tumor), a un tipo de célula (CD4+/CD8+) o a un punto temporal (pretratamiento/postratamiento) (Fig. S3a). Realizamos el mismo análisis utilizando la distancia de Levenshtein por pares, y los resultados siguieron siendo coherentes (Fig. S3b). Los mismos resultados se mantienen cuando utilizamos incrustaciones de redes neuronales para representar las secuencias de TCR, como TCR-BERT o nuestro autoencoder convolucional (Fig. 3c, Fig. S3c). A nivel del gen, solo los genes que codifican las cadenas de TCR muestran correlaciones positivas débiles (Fig. S3d).

Además, investigamos otro conjunto de datos que contiene datos de secuenciación de ARN de una sola célula y secuenciación de TCR emparejados de 20 pacientes con cáncer colorrectal (CRC) y 2 pacientes con cáncer duodenal^[27]. Este conjunto de datos incluye muestras de tejido tumoral, tejido normal adyacente y sangre periférica recogidas antes y después del tratamiento con anti-PD-1. Para mantener la coherencia con nuestro análisis principal, nos centramos en las muestras de tumor y sangre periférica y restringimos el análisis a las células T con secuenciación de ARN y TCR emparejadas (Fig. S4a, Métodos). Esto dio como resultado 201.048 células T con perfiles de expresión de ARN medidos en 36.027 genes, cada una emparejada con una secuencia de cadena β de CDR3 de TCR correspondiente. Como en el conjunto de datos de CHCC, las células T que comparten clones de TCR idénticos estaban más cerca en el espacio de expresión de ARN que los pares de células aleatorias (Fig. S4b). Además, las correlaciones de distancia ARN-TCR siguieron siendo débiles incluso cuando se utilizaron incrustaciones de TCR-BERT o autoencoder CNN (Fig. S4c), lo que refuerza la aplicabilidad más amplia de nuestro marco y modelo de pseudoclones.

Basándonos en estos hallazgos, evaluamos a continuación TRIM en un contexto de enfermedad más amplio utilizando un gran conjunto de datos de cáncer con datos de secuenciación de ARN y TCR emparejados^[11]. Filtramos el conjunto de datos para los pacientes que tenían tanto una muestra de tejido normal adyacente como una muestra de tumor, cada una con al menos 1000 células, y para un tipo de cáncer que tenía al menos cinco pacientes (Fig. S5a). El conjunto de datos filtrado final incluye 19 pacientes de 3 tipos de cáncer (cáncer de esófago, cáncer de tiroides, cáncer de endometrio) y 87.304 células T. A pesar de algunas variaciones entre los tipos de cáncer, observamos constantemente que los clones de TCR idénticos están más cerca en el espacio de ARN que las células aleatorias (Fig. S5b, c), de forma similar a los conjuntos de datos de CHCC y CRC. Además, la falta de correlación entre las distancias de ARN y TCR persiste (Fig. S5d), lo que apoya la aplicabilidad más amplia de nuestro concepto y modelo de pseudoclones. En resumen, en múltiples entornos confirmamos que existe poca relación entre la similitud de la secuencia de TCR y la expresión de ARN: los clones con secuencias de TCR muy diferentes pueden presentar perfiles de expresión de ARN similares.

A continuación, probamos si esto también es cierto para los clones grandes, donde esperamos tener células en diferentes estados de activación. Hipotetizamos que incluso dentro de las regiones del espacio de expresión de ARN donde hay células de alta clonalidad, cualquier clon individual tendrá células cuya expresión de ARN esté más cerca entre sí que de las de otros clones. Para evaluar esta hipótesis, utilizamos el conjunto de datos de CHCC. Como se mencionó anteriormente, las células con recuentos de clones altos (mayores de 10) están enriquecidas en el espacio de ARN correspondiente a las células T CD8+, mientras que las células con recuentos de clones no mayores de 10 se distribuyen por todo el espacio de ARN (Fig. 4a). Calculamos la distancia media por pares de ARN entre las células del mismo clon de TCR (distancia intra-clon) y la comparamos con la distancia media de ARN entre las células de diferentes clones (distancia inter-clon) para los clones de células T CD8+ de alto recuento en la sangre (Fig. 4b) y en el tumor (Fig. 4c, d). Para las células T infiltrantes tumorales, distinguimos entre células citotóxicas y no citotóxicas basándonos en la firma génica informada en la ref. ^[19] (promediando la expresión de PRF1, GZMA, GZMB, GZMH, NKG7, GNLY). En todos los casos, las distancias intra-clon fueron significativamente menores que las distancias inter-clon en promedio en el espacio de ARN (prueba t unilateral, p < 0,01; Fig. 4b-d). Para garantizar que esta relación no se debiera únicamente a que las células provengan de la misma fuente de tejido (sangre o tumor) o compartan el mismo estado, repetimos el análisis dentro de los grupos de células individuales. El mismo patrón se mantuvo en estos subconjuntos más detallados (Fig. S6a). Por último, ampliamos el análisis al pequeño número de células T CD4+ de alta clonalidad en el tumor y observamos la misma tendencia de mayor distancia inter-clon en comparación con la distancia intra-clon de ARN (Fig. S6b).

Por lo tanto, nuestro análisis demuestra que las células T de alta clonalidad, especialmente las células CD8+, exhiben una fuerte relación entre la clonalidad del TCR y la expresión de ARN. Las células del mismo clon están constantemente más cerca unas de otras en el espacio de ARN que las células de diferentes clones, incluso dentro de grupos detallados.

TRIM modela con precisión el espacio de distribución de ARN y la clonalidad del TCR

Para evaluar la capacidad de nuestro modelo para integrar las modalidades de ARN y TCR, comparamos los datos generados fuera de la muestra con los datos de referencia. Primero, evaluamos TRIM en el conjunto de datos de CHCC de cabeza y cuello. La salida del modelo en el espacio de ARN generalmente coincide con los datos de referencia; esto se demuestra para cada distribución condicional, para los pacientes excluidos (Figura S7a). Es importante destacar que, como se muestra en la Figura 5a, el UMAP del espacio de ARN para las células de referencia, coloreado por el recuento real de clones de TCR, coincide visualmente con el UMAP de las células de ARN generadas y su recuento de pseudoclones correspondiente. Para evaluar si TRIM captura la estructura del espacio de ARN, realizamos un análisis de agrupamiento. Específicamente, agrupamos las células en el espacio de expresión de ARN utilizando k-means (con k = 10, 25 y 50), calculamos el recuento medio de clones en cada grupo tanto para los datos reales como para los generados, y luego calculamos la correlación entre estos valores. En otras palabras, dividimos el espacio de ARN en grupos distintos y evaluamos si TRIM reproduce los patrones correspondientes en el espacio de TCR. En todos los niveles de resolución de agrupamiento, las correlaciones fueron consistentemente altas (r = 0,92, 0,95 y 0,92 para k = 10, 25 y 50, respectivamente; Figura S7b). Esto demuestra que nuestro modelo captura con precisión el espacio conjunto a nivel de distribución, y que las discretizaciones de pseudoclones son comparables con los clones reales inherentemente discretos. Si bien estas son características globales de todas las células T del conjunto de datos, a continuación, evaluamos si nuestro modelo también modela con precisión segmentos más granulares de los datos, como las distribuciones marginales a lo largo de cada condición.

A continuación, examinamos la diversidad de clonotipos capturada por TRIM para cada condición, cuantificada mediante la entropía de Shannon (Métodos). La diversidad clonal real y la diversidad de pseudoclones predicha correspondiente estuvieron fuertemente correlacionadas (r = 0,74; Figura 5b). Para evaluar aún más esta relación, también medimos la diversidad como la proporción de TCR únicos en cada muestra, calculada dividiendo el número de TCR distintos por el número total de células. Esta métrica alternativa arrojó una correlación igualmente alta entre los valores reales y predichos (r = 0,80; Figura S7c). En conjunto, estos resultados indican que TRIM captura eficazmente las características clave del repertorio de TCR dentro de su espacio de incrustación aprendido, tanto en las fuentes de tejido (sangre y tumor) como en los puntos de tiempo del tratamiento (pretratamiento y postratamiento).

Luego, evaluamos TRIM en el conjunto de datos de CRC para evaluar su capacidad de generalización. Las salidas del modelo en el espacio de ARN coincidieron estrechamente con los datos de referencia, y la correlación entre la diversidad de clonotipos reales y predichos fue alta (r = 0,67; Figura S8a, b). Un análisis de agrupamiento más detallado también demostró un fuerte acuerdo entre los datos reales y generados, con altas correlaciones en los recuentos medios de clones de TCR en los diferentes niveles de resolución de agrupamiento (r = 0,94, 0,85 y 0,86 para k = 10, 25 y 50, respectivamente; Figura S8c).

Para examinar la robustez de TRIM en diferentes poblaciones de pacientes, investigamos factores como la edad y el sexo utilizando el conjunto de datos de CRC. También exploramos cómo la variación del número de pacientes en los datos de entrenamiento afectó el rendimiento del modelo. Específicamente, realizamos un análisis de exclusión de grupos, donde se excluyeron grupos completos durante el entrenamiento (por ejemplo, excluyendo a todas las pacientes, o a los pacientes en el percentil superior o inferior del 15 para la edad o el recuento de muestras), y utilizamos estos grupos para las pruebas (Figura S9a). En estas diferentes condiciones, TRIM capturó consistentemente el espacio de ARN-TCR conjunto a nivel de distribución, como lo indican las altas correlaciones en los diferentes niveles de resolución de agrupamiento (k = 10, 25 y 50; Figura S9b-f). Observamos que un agrupamiento de mayor resolución dificulta la predicción de los recuentos medios de clones de TCR, y esta dificultad se amplifica aún más al evaluar muestras con pocas células, ya que los grupos más pequeños contienen información muy limitada para una estimación precisa. Sin embargo, TRIM mantuvo un fuerte rendimiento en la predicción de la diversidad de clonotipos, con altas correlaciones entre los valores reales y predichos (Figura S9b-f, paneles de la derecha). Observamos una reducción de la capacidad de generalización cuando se excluyeron los pacientes de mayor edad del entrenamiento (Figura S9d). Hipotetizamos que esta dificultad en la generalización puede deberse a las propiedades distintas de esta subpoblación de pacientes, que muestra una reducción de la diversidad clonal de las células T. No obstante, los resultados de todos los análisis de subgrupos fueron consistentes, lo que indica que TRIM es robusto a las variaciones demográficas y del tamaño de la muestra.

También evaluamos nuestro modelo en el conjunto de datos de cáncer generalizado ^[11] (Figura S10a). Para evaluar cuantitativamente las predicciones de TRIM, agrupamos las células en el espacio de ARN utilizando k-means (k = 20) y luego calculamos la diversidad clonal real y predicha y el recuento medio de clones para ese grupo (Figura S10b). Encontramos un alto grado de correlación en ambos casos, r = 0,97 y r = 0,97, respectivamente. Además, las proyecciones UMAP de los datos reales y los datos predichos de nuestro modelo mostraron que el modelo había aprendido a generar, para los pacientes excluidos, perfiles de expresión de ARN precisos basados en la fuente (tejido normal/tumor) y el tipo de cáncer (cáncer de esófago o ESCA, cáncer de tiroides o THCA, carcinoma de cuerpo uterino o UCEC) (Figura S10c, gráficos de la izquierda). Además, el modelo también aprendió a igualar el recuento de pseudoclones en cada región del espacio de ARN con el recuento de clones de referencia correspondiente para esa región (Figura S10c, gráficos de la derecha). En conjunto, estos resultados indican que el modelo logra una alta precisión tanto a nivel del conjunto de datos como a nivel del grupo.

TRIM predice con precisión la expansión clonal y el estado de activación de las células T después del tratamiento

A continuación, consideramos la resolución más granular: la predicción de la expansión a nivel de células individuales. En el conjunto de datos de CHCC de cabeza y cuello, las muestras de sangre previas al tratamiento incluyeron 15.147 clones de TCR únicos. De estos, 725 persistieron en la sangre postratamiento, y 185 se expandieron postratamiento. En particular, 66 de estos 185 clones expandidos (36%) también se detectaron en el tumor postratamiento. Y 64 de estos 66 TCR se expandieron en el tumor postratamiento. Esta observación respaldó los hallazgos previos de una conexión entre la dinámica de las células T periféricas (sangre) y las asociadas al tumor ^[19], lo que nos motivó a modelar la expansión de las células T utilizando solo los datos derivados de la sangre. Debido al desafío de predecir directamente las secuencias de TCR, como hemos demostrado anteriormente, que las secuencias de ARN y TCR de una sola célula están relacionadas solo a través del recuento de clones, modelamos la expansión a nivel de clon y restringimos el análisis a los TCR presentes en ambas condiciones. Dado un TCR específico y su recuento de clones en la sangre pretratamiento, nuestro objetivo era predecir si su recuento de clones aumentaría o disminuiría postratamiento (normalizado por el número total de células en cada condición). Para hacer esto, mapeamos el TCR al espacio latente y obtuvimos el vector latente correspondiente z. Luego, reemplazamos la incrustación condicional aprendida para "pretratamiento" y, en cambio, utilizamos la incrustación aprendida para "postratamiento", y decodificamos como de costumbre. Comparamos el recuento de pseudoclones resultante con el recuento de pseudoclones obtenido decodificando los vectores latentes originales, y si el primero era mayor, predijimos la expansión. Definimos la expansión en los datos de referencia como un TCR que se observa en ambos puntos de tiempo y cuyo recuento en la distribución de postratamiento es mayor que en la distribución de pretratamiento, normalizado por el número total de células en cada condición (consulte Métodos para obtener más detalles). Como se muestra en la Figura 5c, nuestro modelo logró un área bajo la curva de características operativas del receptor (AUC) de 0,92 para esta tarea de predicción. En contraste, los métodos de referencia, como un clasificador de vecinos más cercanos (KNN) o un perceptrón multicapa (MLP), tuvieron un rendimiento similar a la línea de referencia aleatoria (línea diagonal) al predecir una etiqueta binaria de expansión, cuando se les proporcionó como entrada la expresión de ARN, las incrustaciones del codificador automático de TCR o la concatenación de ambos.

De manera similar, evaluamos TRIM en el conjunto de datos de CRC. En la evaluación de exclusión de un elemento, TRIM logró un AUC de 0,83 para predecir la expansión clonal postratamiento, mientras que los modelos de referencia (KNN y MLP) tuvieron un rendimiento inferior o cercano al aleatorio (incluso después de una optimización exhaustiva de hiperparámetros), independientemente de si utilizaron la expresión de ARN, las incrustaciones de TCR o su combinación como entrada (Figura S11a-d, Métodos). Para evaluar aún más la robustez, probamos TRIM en entornos de exclusión de grupos, donde los pacientes de prueba tenían mediciones de postratamiento. TRIM mantuvo un rendimiento competitivo, logrando AUC de 0,68, 0,81 y 0,78 para las evaluaciones de exclusión de mujeres, exclusión de jóvenes y exclusión de más células, respectivamente, mientras que los modelos de referencia nuevamente tuvieron un rendimiento cercano al aleatorio (Figura S12).

Luego, evaluamos la capacidad de TRIM para predecir la expansión clonal a nivel de célula individual en el conjunto de datos de cáncer generalizado (Figura S10d). Predijimos si un TCR en particular que se encuentra en el tejido normal adyacente tenía un recuento de clones más alto en el tumor. Realizamos esta tarea de predicción de la misma manera generativa condicional que en el conjunto de datos de CHCC de cabeza y cuello, restringiéndonos solo a los TCR que aparecen en ambas condiciones.

Como antes, los clasificadores KNN y las redes neuronales de referencia que recibieron la concatenación de las representaciones de ARN y TCR tuvieron un rendimiento similar a la línea de referencia aleatoria, mientras que nuestro modelo superó claramente a estas líneas de referencia. Sin embargo, el AUC es menor que en el conjunto de datos de CHCC de cabeza y cuello (0,57 frente a 0,92), lo que puede explicarse por el hecho de que las poblaciones de células T en la sangre periférica circundante contienen más información sobre las poblaciones de células T en el tumor que las que se encuentran en el tejido normal adyacente. Hipotetizamos que la capacidad del modelo para predecir la expansión de las células T en el tumor a partir del tejido normal adyacente mejor que al azar se debe a la señalización y la comunicación entre estas poblaciones de células T, aunque en menor medida que con las células T circulantes en la sangre.

Para validar la utilidad de nuestro modelo para los análisis biológicos posteriores, comparamos sus predicciones basadas en los datos de sangre pretratamiento con los análisis de los fenotipos de tumores y células T postratamiento informados en la referencia ^[19] en el conjunto de datos de CHCC de cabeza y cuello (Figura 5d, e, Figura S13). En ^[19], los autores informaron sobre la expansión de las células T CD8 proliferativas CD38+/HLA-DR+ en la sangre durante el curso del tratamiento y encontraron que la frecuencia de las células T CD8 activadas KLRG1-negativas PD-1+ en la sangre pretratamiento y en el tratamiento estaba fuertemente asociada con la respuesta patológica intratumoral. Primero, calculamos el porcentaje de células T proliferativas CD38+/HLA-DR+ para cada condición en función de nuestros datos predichos y encontramos que esta población aumenta en las condiciones de postratamiento (Figura S13). A nivel de paciente y punto de tiempo, las predicciones de nuestro modelo están bien correlacionadas con los datos de referencia para esta estadística (r = 0,89). El parámetro de pendiente de una regresión de los valores reales de esta métrica con los valores predichos inferidos a partir de la salida de nuestro modelo también se muestra, lo que da como resultado β1 = 0,82. A continuación, modelamos la expresión de KLRG1 en las células CD8 PD-1+ y encontramos que disminuye más para los pacientes que tienen una respuesta patológica que para aquellos que no la tienen, con una correlación por paciente de 0,82 (Figura 5d). Por último, aislamos las células CD8+ predichas, calculamos la proporción que son CD38+/HLA-DR+ y encontramos que esto también está fuertemente correlacionado con la respuesta patológica (Figura 5e). Las predicciones de nuestro modelo están muy correlacionadas con los datos reales por paciente (r = 0,90). En general, estos análisis muestran que TRIM genera predicciones precisas y significativas de la expansión clonal.

TRIM identifica genes asociados con la expansión clonal

Para analizar más a fondo el aspecto de la predicción de expansión de nuestro modelo, investigamos las características que utiliza para realizar esta tarea. Específicamente, realizamos un análisis de expresión génica diferencial al comparar todos los pseudoclones predichos para expandirse por nuestro modelo con aquellos que no se predijo que se expandieran en las células T de la sangre. Luego, clasificamos los genes según su cambio logarítmico en la expresión, sujeto a la restricción de que se expresen en al menos el 10% de las células. Repetimos este proceso para los datos reales, basándonos en la expansión observada y su expresión génica medida en cada condición. Luego, comparamos las dos listas clasificadas considerando la intersección de sus primeros k elementos, para valores crecientes de k. La superposición de un modelo aleatorio estaría en la línea diagonal, obteniendo solo una intersección de genes en promedio. Primero, evaluamos TRIM en el conjunto de datos de CHCC de cabeza y cuello. El AUC de nuestro modelo con respecto a esta línea de base es de 0,92 (Fig. 6a). Repetimos este proceso por separado para los subconjuntos de células T CD4+ y CD8+, lo que resultó en AUC de 0,80 y 0,77, respectivamente (Figura suplementaria 14a). También evaluamos este aspecto de nuestro modelo tomando una firma génica de expansión y una firma génica de no expansión que se identificaron previamente en la referencia ^[19]. Cuando restringimos nuestro análisis a estos genes y los clasificamos por el cambio logarítmico en la expresión, los genes de la firma de expansión tenían más probabilidades de estar sobreexpresados en los clones que el modelo predijo que se expandirían (Fig. 6a). La cuantificación utilizando una prueba de rango de Wilcoxon ^[45] para la hipótesis de que el modelo clasifica los genes de la firma de expansión más alto que los genes de la firma de no expansión resultó en un valor p de 0,064. Dado que las firmas se derivaron de tejido tumoral, mientras que nuestro análisis modela la expansión en la sangre, esperamos cierta atenuación; no obstante, el resultado muestra un enriquecimiento direccionalmente consistente y se alinea con la evidencia previa de que los programas de expansión asociados al tumor son detectables en la sangre. Esto demuestra que nuestro modelo predice con precisión no solo qué clones se expandirán, sino también qué genes se expresarán de manera diferente durante esa expansión.

Llevando este análisis un paso más allá, investigamos el modelo para identificar los genes que se predice que se expresarán de manera diferente (DE) desde el pretratamiento hasta el postratamiento por los pseudoclones predichos para expandirse en la sangre (Figura suplementaria 14b). El análisis de enriquecimiento de conjuntos de genes de los genes DE reveló un enriquecimiento para las vías del ciclo celular y la proliferación celular, que se espera que estén involucradas en la expansión clonal (Figura suplementaria 14c). Estos resultados respaldan la precisión de las predicciones de expresión génica de nuestros pseudoclones.

Luego, evaluamos la capacidad de TRIM para identificar genes expresados de manera diferente en el conjunto de datos de CRC. En la evaluación de eliminación de uno, TRIM logró un AUC de 0,85 para todas las células T, 0,82 para las células CD4 y 0,79 para las células CD8 (Figura suplementaria 15a). En diferentes análisis de subgrupos, TRIM mantuvo un fuerte rendimiento, con AUC para todas las células T de 0,83, 0,90 y 0,81 para la eliminación de mujeres, la eliminación de jóvenes y la eliminación de más células, respectivamente; para las células CD4, 0,83, 0,83 y 0,81; y para las células CD8, 0,72, 0,84 y 0,72 (Figura suplementaria 15b–d). Estos resultados validan aún más la capacidad de TRIM para generalizarse a diferentes tipos de cáncer, así como en diversas condiciones de los pacientes.

Para identificar los genes que el modelo utiliza para predecir la expansión clonal de las células T CD8, realizamos un análisis de saliencia. Esto implicó la creación de un clasificador para la expansión clonal basado en el espacio latente de la red y el examen de los gradientes en relación con cada gen de entrada. A través de este análisis de saliencia, generamos una clasificación de genes basada en su poder predictivo para la expansión, que va desde aquellos que son fuertemente predictivos negativos hasta aquellos que son fuertemente predictivos positivos. Comparamos esta clasificación (gradiente) con los resultados de los análisis DE de rutina para evaluar su relevancia. Específicamente, examinamos tres escenarios DE: (1) grupos expandidos frente a no expandidos en el postratamiento, (2) grupos que se expanden frente a no expandidos en el pretratamiento y (3) un caso simulado sin datos de postratamiento, comparando células con tamaños de clon grandes frente a tamaños de clon pequeños. Para cada método, identificamos los genes importantes en función de los valores p corregidos (para los análisis DE) o las puntuaciones de clasificación (en nuestro método). Los análisis de correlación de Pearson mostraron que los genes identificados como importantes para que nuestro modelo prediga la expansión clonal se alinean más estrechamente con la comparación de expansión en el pretratamiento frente a no expansión (0,52) que con los otros dos enfoques DE (-0,07 y 0,22) (Fig. 6b). De la Fig. 6c, observamos que los genes de gradiente tienen una superposición limitada con los genes identificados a través de los cambios logarítmicos en la expresión en las comparaciones de postratamiento y tamaño de clon, y un número relativamente mayor de genes compartidos con la lista DE de pretratamiento.

Para evaluar aún más la relevancia biológica de los genes identificados por cada enfoque, realizamos un análisis de enriquecimiento de términos de la ontología génica (GO) utilizando el paquete gseapy ^[46] (Figura suplementaria 16, Métodos). Descubrimos que los conjuntos de genes identificados como enriquecidos entre los genes basados en el gradiente coincidieron en gran medida con los conjuntos de genes identificados utilizando los análisis DE más estándar (Figura suplementaria 16). Luego, consultamos a ChatGPT para compilar una lista de genes asociados con la expansión de las células T y analizamos las listas resultantes en los diferentes métodos (Fig. 6d). Este análisis se inspiró en los hallazgos de la referencia ^[47], que demostró el potencial de los modelos de lenguaje grandes para anotar conjuntos de genes. Nuestro método ("Gradiente"), junto con los métodos supervisados que utilizaron las etiquetas de postratamiento ("Expandido frente a no expandido" y "Expansión frente a no expansión"), identificó con éxito los genes involucrados en el ciclo celular (AURKA, CENPU, DUT), la función efectora de las células T y la citotoxicidad (PRF1, GZMA, GZMH, GNLY, NKG7, CST7), la activación y la coestimulación de las células T (CTLA4, TNFRSF9), las vías metabólicas y mitocondriales (ATP5PF, VDAC2, DHTKD1, IBA57) y la migración de las células T (CCR7, CCL5, IL6R, FOXO1). En particular, las células T CD38+ CD8 que se expandieron en la sangre durante el tratamiento se ha demostrado que exhiben firmas de activación y proliferación; aquí, la identificación de una firma de migración de células T es nueva y sugiere que estas células pueden desplazarse hacia el tumor y los ganglios linfáticos, lo que es consistente con que esta población tenga TCR que se expandieron en los tumores durante la terapia. Es importante destacar que, si bien los análisis de expresión diferencial (DE) supervisados sirven como límites superiores teóricos, porque se basan directamente en las etiquetas de postratamiento, nuestro método pudo recuperar señales biológicas similares sin requerir mediciones de postratamiento en el momento de la inferencia. Aunque TRIM requiere acceso a datos emparejados de pretratamiento y postratamiento durante el entrenamiento, una vez entrenado, se puede aplicar a nuevas muestras de solo pretratamiento para identificar programas génicos asociados con la expansión de las células T. Esta capacidad proporciona tanto la interpretabilidad del modelo, al demostrar que sus predicciones se basan en señales biológicas significativas, como la relevancia clínica, ya que proporciona estas perspectivas predictivas en un entorno clínico común donde las muestras de postratamiento no están disponibles. Esto permite que TRIM proporcione información significativa en entornos donde las mediciones de postratamiento no están disponibles, lo que cierra la brecha entre los análisis supervisados y la aplicabilidad en el mundo real. En contraste, los genes expresados de manera diferente identificados al comparar clones grandes frente a clones pequeños en el pretratamiento se enriquecieron para los programas citotóxicos, efectores e inmunorreguladores, pero no mostraron un enriquecimiento significativo para los programas relacionados con la expansión (evaluado mediante pruebas hipergeométricas por conjunto de genes con corrección FDR; el enriquecimiento se muestra como valores p ajustados por FDR en la Fig. 6d). Esto sugiere que TRIM proporciona una visión más informativa y biológicamente fundamentada de la expansión de las células T en condiciones del mundo real.

TRIM utiliza un autoencoder variacional condicional para predecir el efecto del tratamiento en el tumor a partir de la sangre en el pretratamiento

Para integrar la información de las secuencias de ARN y TCR en las células T y realizar predicciones sobre los estados de las células T en sitios y puntos de tiempo no vistos, proponemos TRIM, un modelo generativo basado en un VAE condicional (Fig. 1b, Métodos: TRIM). TRIM se entrena con datos de secuenciación de ARN y TCR de células individuales de células T en la sangre circulante antes del tratamiento, y luego predice el estado de la población de células T a nivel de célula individual tanto en los tumores de pretratamiento como de postratamiento, así como en las muestras de sangre de postratamiento (Fig. 1b). TRIM aprende una representación conjunta de la expresión génica y la clonalidad utilizando dos codificadores separados para mapear el perfil de ARN y la secuencia de TCR de una célula en un espacio latente conjunto, condicionado a la fuente (sangre o tumor), el punto de tiempo (pretratamiento o postratamiento) y el paciente. Para garantizar que el espacio latente conjunto capture información de ambos, el perfil de ARN de la célula y su secuencia de TCR, TRIM utiliza dos decodificadores separados para reconstruir las dos modalidades (Fig. 2a).

Cuando están presentes múltiples modalidades emparejadas, los autoencoders tradicionales minimizan el error al traducir entre modalidades utilizando el codificador de una y el decodificador de la otra. Esta estrategia ha sido efectiva para integrar varias modalidades de células individuales, como RNA-seq, ATAC-seq e imágenes de cromatina ^[35], ^[38], ^[39]. Se aplica una pérdida de reconstrucción para alentar a que la salida se acerque a los puntos originales. Para el espacio de ARN, esta pérdida se calcula utilizando la norma L2 (Métodos). Las representaciones de TCR utilizadas en nuestro modelo se derivan de incrustaciones aprendidas a través de un autoencoder convolucional, que entrenamos como un paso de preprocesamiento. Específicamente, cada TCR está representado por la secuencia de aminoácidos CDR3 de su cadena β, codificada como una secuencia con relleno nulo donde cada aminoácido está representado por sus factores de Atchley, descriptores numéricos informados bioquímicamente que capturan propiedades como la polaridad, la carga y la hidrofobicidad ^[40] (Métodos). Extraemos las incrustaciones de la capa latente media del autoencoder convolucional entrenado para que sirvan como entrada para el modelo TRIM. A pesar de centrarse en la secuencia CDR3 altamente variable de las cadenas β de los TCR, el modelo aprende implícitamente los segmentos génicos V y J contribuyentes, lo que da como resultado un espacio de representación aprendido organizado de manera biológicamente significativa (Figura suplementaria 1). Nos centramos en la cadena β por razones prácticas y biológicas. Primero, en la secuenciación V(D)J de células individuales, el emparejamiento α-β a menudo es incompleto. En las muestras de sangre perfiladas en la referencia ^[19], el 87% de las células tienen cadenas emparejadas, con β recuperada con más frecuencia que α (95% frente a 89%). En general, los grandes recursos públicos de TCR son predominantemente solo β, con el 94% de las secuencias que contienen solo la cadena β ^[41]. En segundo lugar, los estudios anteriores que incluyen ambas cadenas informan solo ganancias marginales en el rendimiento, probablemente porque la secuencia CDR3β es el determinante dominante de la especificidad péptido-MHC ^[42], ^[43]. Dadas estas limitaciones experimentales y la disponibilidad de datos, utilizamos la cadena β para maximizar la cobertura de la muestra y, al mismo tiempo, conservar la mayor parte de la señal relevante.

También observamos que los perfiles de RNA-seq contienen información limitada sobre la secuencia de TCR de una célula más allá de su clonalidad (como discutiremos con más detalle más adelante), una traducción directa entre TCR y RNA no es posible. Para abordar esta limitación, modelamos la clonalidad de las células T en lugar de las secuencias de TCR. Para preservar la estructura local en el espacio de TCR sin penalizar las variaciones no relacionadas, introducimos la pérdida de clonalidad contrastiva. Esta pérdida alienta a que las células T del mismo clon estén cerca en el espacio latente, al tiempo que garantiza que las diferentes células T se mantengan separadas por un margen (Métodos, pérdida de clonalidad contrastiva). Este enfoque cambia el enfoque de la representación exacta de TCR a la preservación de la clonalidad.

Además, nuestro modelo incorpora las condiciones del paciente, la fuente y el punto temporal como entradas, representadas por vectores de longitud fija que se aprenden conjuntamente con el resto del modelo a través de una red de incrustación separada, F (Fig. 2a). En contraste con el enfoque estándar, que aprende vectores condicionales indexando una matriz de incrustación inicializada aleatoriamente (con los parámetros de esta matriz aprendidos durante el entrenamiento), nuestro método está diseñado para aprender las representaciones vectoriales de las condiciones utilizando los perfiles de ARN y TCR asociados con cada condición. Específicamente, la red F procesa un subconjunto de células de entrenamiento, Semb, y las mapea en un espacio latente de demb dimensiones. Para calcular el vector de incrustación para una condición específica, TRIM promedia las representaciones de las células que corresponden a esa condición. Por ejemplo, el vector del punto temporal pretratamiento se calcula promediando las representaciones de las células pretratamiento en Semb (Métodos, Aprendizaje de incrustaciones de condiciones). Para las incrustaciones de pacientes, nos centramos específicamente en las células sanguíneas pretratamiento, ya que estos son los datos disponibles para los pacientes no vistos. Las incrustaciones se aprenden a partir de los datos de ARN y TCR, utilizando dos redes, FRNA y FTCR, cuyas salidas se promedian para formar la incrustación final. Es importante destacar que nuestro modelo puede aprovechar el hecho de que la distribución de ARN o TCR de un paciente de una fuente/punto temporal es muy informativa sobre sus distribuciones en otros entornos. Demostramos que este enfoque es más eficaz para generar distribuciones precisas de ARN y TCR en diferentes condiciones y, lo que es crucial, permite una extrapolación precisa a nuevos pacientes no vistos. Durante la evaluación del modelo, discretizamos las salidas continuas del decodificador para estimar los clones de TCR, asignando pseudoclonos utilizando un algoritmo voraz que agrupa los puntos dentro de un margen en el mismo clon (Fig. 2b, Métodos, Asignación de pseudoclonos). Para maximizar la utilidad de los datos y probar su generalizabilidad, entrenamos con todos los datos excepto los de un paciente y probamos el modelo con el paciente excluido. Repetimos este procedimiento en un proceso de validación cruzada para cada paciente, con el fin de obtener el conjunto completo de predicciones. Esto permite comparar la salida con la verdad fundamental en todo el conjunto de datos, pero todo proviene de predicciones fuera de la muestra.

Las secuencias de ARN y TCR de una sola célula están relacionadas solo a través del recuento de clones

Para explorar la relación entre TCR y ARN, analizamos datos emparejados de scRNAseq y scTCRseq de 27 pacientes con CHCC^[19]. El conjunto de datos incluye mediciones tanto de sangre periférica como de células T infiltrantes tumorales tomadas antes y después de la terapia con inhibidores de puntos de control inmunitario. Consta de 73.961 células T con expresión de ARN medida en 29.750 genes, emparejadas con secuencias de la cadena β de CDR3 de TCR. Si bien cada perfil de expresión de ARN es único, muchas células son clones con secuencias de TCR idénticas, y solo el 69,8% de las células T tienen TCR únicos. Si bien estudios previos sugieren que los clones de células T idénticos comparten expresiones de ARN similares^[32], aún quedan muchas preguntas sobre la relación entre el espacio de expresión de ARN y el espacio de TCR.

Para ello, primero visualizamos los datos de scRNAseq utilizando UMAP^[44] (Fig. 3a). Se puede observar un fuerte solapamiento en las distribuciones de cada paciente (Fig. 3a^[44], esquina inferior izquierda), con cuatro grupos distintos: las células T de la sangre están separadas de las del tumor, y las células T CD4+ son distintas de las células T CD8+ en ambas fuentes. Si bien los tipos de células están claramente separados en el espacio transcriptómico, el efecto del tratamiento es más sutil. No se observa un cambio discernible en las células sanguíneas, pero se produce un ligero cambio en las células tumorales, particularmente en las células CD8+. Al asignar a cada célula su recuento de clones (que corresponde al número de células que expresan el mismo TCR), se puede observar una mayor clonalidad en las regiones del espacio de expresión de ARN que corresponden a partes de los grupos de células T CD8+. Los tamaños de los clones más grandes son más frecuentes en las células CD8+ y en el tumor en comparación con las células CD4+ y las de la sangre (Figura suplementaria 2).

Para evaluar la hipótesis de que los clones están más cerca en la expresión de ARN, comparamos las células por distancia por pares en el espacio de ARN frente a la distancia de Hamming por pares en el espacio de TCR (Fig. 3b). La distribución de las distancias de ARN entre los clones de TCR idénticos (distancia de Hamming de cero) es sustancialmente menor que las distancias de ARN para las células que están más separadas en el espacio de TCR. Curiosamente, una vez que superamos los clones de TCR idénticos, existe muy poca relación entre la secuencia de TCR y la expresión de ARN de dos células. Para encontrar alguna relación, debemos filtrar los clones con un recuento mayor que uno, e incluso entonces, la diferencia entre la distancia media de ARN entre los clones separados por 0 de Hamming y los separados por 1 de Hamming es 4 veces mayor que la diferencia entre las medias de las células separadas por 1 de Hamming y 10 de Hamming. En otras palabras, dos células cuyas secuencias de TCR difieren en solo una posición están tan lejos en el espacio de ARN como dos células con secuencias de TCR completamente diferentes. Observamos que esto es cierto no solo a nivel del conjunto de datos, sino también si se restringe solo a un paciente, a una fuente (sangre/tumor), a un tipo de célula (CD4+/CD8+) o a un punto temporal (pretratamiento/postratamiento) (Figura suplementaria 3a). Realizamos el mismo análisis utilizando la distancia de Levenshtein por pares, y los resultados siguieron siendo consistentes (Figura suplementaria 3b). Los mismos hallazgos se mantienen cuando utilizamos incrustaciones de redes neuronales para representar las secuencias de TCR, como TCR-BERT o nuestro autoencoder convolucional (Fig. 3c, Figura suplementaria 3c). A nivel del gen, solo los genes que codifican las cadenas de TCR muestran correlaciones positivas débiles (Figura suplementaria 3d).

Además, investigamos otro conjunto de datos que contiene datos emparejados de ARN-seq y TCR-seq de una sola célula de 20 pacientes con cáncer colorrectal (CCR) y 2 pacientes con cáncer duodenal^[27]. Este conjunto de datos incluye muestras de tejido tumoral, tejido normal adyacente y sangre periférica recolectadas antes y después del tratamiento anti-PD-1. Para mantener la coherencia con nuestro análisis principal, nos centramos en las muestras de tumor y sangre periférica y restringimos el análisis a las células T con datos emparejados de ARN-seq y TCR-seq (Figura suplementaria 4a, Métodos). Esto resultó en 201.048 células T con perfiles de expresión de ARN medidos en 36.027 genes, cada una emparejada con una secuencia de cadena β de CDR3 de TCR correspondiente. Como en el conjunto de datos de CHCC, las células T que comparten clones de TCR idénticos estaban más cerca en el espacio de expresión de ARN que los pares de células aleatorias (Figura suplementaria 4b). Además, las correlaciones de distancia ARN-TCR siguieron siendo débiles incluso cuando se utilizaron incrustaciones de TCR-BERT o autoencoder CNN (Figura suplementaria 4c), lo que refuerza la aplicabilidad más amplia de nuestro marco y modelo de pseudoclones.

Basándonos en estos hallazgos, evaluamos a continuación TRIM en un contexto de enfermedad más amplio utilizando un gran conjunto de datos de cáncer con datos de secuenciación de ARN y TCR emparejados^[11]. Filtramos el conjunto de datos para los pacientes que tenían tanto una muestra de tejido normal adyacente como una muestra de tumor, cada una con al menos 1000 células, y para un tipo de cáncer que tenía al menos cinco pacientes (Figura suplementaria 5a). El conjunto de datos filtrado final incluye 19 pacientes de 3 tipos de cáncer (cáncer de esófago, cáncer de tiroides, cáncer de endometrio) y 87.304 células T. A pesar de cierta variación entre los tipos de cáncer, observamos constantemente que los clones de TCR idénticos están más cerca en el espacio de ARN que las células aleatorias (Figura suplementaria 5b, c), de manera similar a los conjuntos de datos de CHCC y CCR. Además, la falta de correlación entre las distancias de ARN y TCR persiste (Figura suplementaria 5d), lo que respalda la aplicabilidad más amplia de nuestro concepto y modelo de pseudoclonos. En resumen, en múltiples entornos confirmamos que existe poca relación entre la similitud de la secuencia de TCR y la expresión de ARN: los clones con secuencias de TCR muy diferentes pueden presentar perfiles de expresión de ARN similares.

A continuación, probamos si esto también es cierto para los clones grandes, donde esperamos tener células en diferentes estados de activación. Hipotetizamos que incluso dentro de las regiones del espacio de expresión de ARN donde hay células de alta clonalidad, cualquier clon individual tendrá células cuya expresión de ARN esté más cerca entre sí que de las de otros clones. Para evaluar esta hipótesis, utilizamos el conjunto de datos de CHCC. Como se mencionó anteriormente, las células con altos recuentos de clones (mayores de 10) están enriquecidas en el espacio de ARN correspondiente a las células T CD8+, mientras que las células con recuentos de clones no mayores de 10 se distribuyen en todo el espacio de ARN (Fig. 4a). Calculamos la distancia media de ARN por pares entre las células del mismo clon de TCR (distancia intra-clon) y la comparamos con la distancia media de ARN entre las células de diferentes clones (distancia inter-clon) para los clones de células T CD8+ de alto recuento en la sangre (Fig. 4b) y en el tumor (Fig. 4c, d). Para las células T CD8+ infiltrantes tumorales, distinguimos entre células citotóxicas y no citotóxicas basándonos en la firma génica informada en la ref. ^[19] (promediando la expresión de PRF1, GZMA, GZMB, GZMH, NKG7, GNLY). En todos los casos, las distancias intra-clon fueron significativamente menores que las distancias inter-clon en promedio en el espacio de ARN (prueba t unilateral, p < 0,01; Fig. 4b–d). Para garantizar que esta relación no se debiera únicamente a que las células provengan de la misma fuente de tejido (sangre o tumor) o compartan el mismo estado, repetimos el análisis dentro de los grupos de estado celular individuales. El mismo patrón se mantuvo en estos subconjuntos más detallados (Figura suplementaria 6a). Por último, ampliamos el análisis al pequeño número de células T CD4+ de alta clonalidad en el tumor y observamos la misma tendencia de mayor distancia inter-clon en comparación con la distancia intra-clon de ARN (Figura suplementaria 6b).

Por lo tanto, nuestro análisis demuestra que las células T de alta clonalidad, especialmente las células CD8+, exhiben una fuerte relación entre la clonalidad de TCR y la expresión de ARN. Las células del mismo clon están consistentemente más cerca entre sí en el espacio de ARN que las células de diferentes clones, incluso dentro de grupos detallados.

TRIM modela con precisión el espacio de distribución de ARN y la clonalidad de TCR

Para evaluar la capacidad de nuestro modelo para integrar las modalidades de ARN y TCR, comparamos los datos generados fuera de la muestra con los datos de la verdad fundamental. Primero evaluamos TRIM en el conjunto de datos de CHCC. La salida del modelo en el espacio de ARN generalmente coincide con la verdad fundamental: esto se muestra en cada distribución condicional, para los pacientes excluidos (Figura suplementaria 7a). Es importante destacar que, como se muestra en la Fig. 5a, el UMAP del espacio de ARN para las células de la verdad fundamental coloreadas por el recuento de clones de TCR reales coincide visualmente con el UMAP de las células de ARN generadas y su recuento de pseudoclon correspondiente. Para evaluar si TRIM captura la estructura del espacio de ARN, realizamos un análisis de agrupamiento. Específicamente, agrupamos las células en el espacio de expresión de ARN utilizando k-means (con k = 10, 25 y 50), calculamos el recuento de clones medio en cada grupo tanto para los datos reales como para los generados y, a continuación, calculamos la correlación entre estos valores. En otras palabras, dividimos el espacio de ARN en grupos distintos y evaluamos si TRIM reproduce los patrones correspondientes en el espacio de TCR. En todas las resoluciones de agrupamiento, las correlaciones fueron consistentemente altas (r = 0,92, 0,95 y 0,92 para k = 10, 25 y 50, respectivamente; Figura suplementaria 7b). Esto demuestra que nuestro modelo captura con precisión el espacio conjunto a nivel de distribución, y que las discretizaciones de pseudoclon son comparables a los clones reales inherentemente discretos. Si bien estas son características globales de todas las células T del conjunto de datos, a continuación, evaluamos si nuestro modelo también modela con precisión segmentos más granulares de los datos, como las distribuciones marginales a lo largo de cada condición.

A continuación, examinamos la diversidad de clonotipos capturada por TRIM para cada condición, cuantificada mediante la entropía de Shannon (Métodos). La verdadera diversidad clonal y la correspondiente diversidad de pseudoclones predicha estuvieron fuertemente correlacionadas (r = 0,74; Fig. 5b). Para evaluar aún más esta relación, también medimos la diversidad como la proporción de TCR únicos en cada muestra, calculada dividiendo el número de TCR distintos por el número total de células. Esta métrica alternativa arrojó una correlación igualmente alta entre los valores reales y predichos (r = 0,80; Fig. S7c). En conjunto, estos resultados indican que TRIM captura eficazmente las características clave del repertorio de TCR dentro de su espacio de incrustación aprendido, tanto en las fuentes de tejido (sangre y tumor) como en los puntos temporales del tratamiento (pretratamiento y postratamiento).

Luego, evaluamos TRIM en el conjunto de datos de CRC para evaluar su capacidad de generalización. Las salidas del modelo en el espacio de ARN coincidieron estrechamente con la verdad fundamental, y la correlación entre la diversidad de clonotipos real y la predicha fue alta (r = 0,67; Fig. S8a, b). Un análisis de agrupamiento más detallado también demostró un fuerte acuerdo entre los datos reales y los generados, con altas correlaciones en los recuentos medianos de clones de TCR en los diferentes niveles de resolución del agrupamiento (r = 0,94, 0,85 y 0,86 para k = 10, 25 y 50, respectivamente; Fig. S8c).

Para examinar la robustez de TRIM en diferentes poblaciones de pacientes, investigamos factores como la edad y el sexo utilizando el conjunto de datos de CRC. También exploramos cómo la variación en el número de pacientes en los datos de entrenamiento afectó el rendimiento del modelo. Específicamente, realizamos un análisis de exclusión de grupos, donde se excluyeron grupos completos durante el entrenamiento (por ejemplo, excluyendo a todas las pacientes, o a los pacientes en el percentil 15 superior o inferior para la edad o el recuento de muestras), y utilizamos estos grupos para las pruebas (Fig. S9a). En estas diferentes condiciones, TRIM capturó consistentemente el espacio de ARN-TCR conjunto a nivel de distribución, como lo indican las altas correlaciones en los diferentes niveles de resolución del agrupamiento (k = 10, 25 y 50; Fig. S9b–f). Observamos que un agrupamiento de mayor resolución dificulta la predicción de los recuentos medianos de clones de TCR, y esta dificultad se amplifica aún más al evaluar muestras con pocas células, ya que los grupos más pequeños contienen información muy limitada para una estimación precisa. Sin embargo, TRIM mantuvo un fuerte rendimiento en la predicción de la diversidad de clonotipos, con altas correlaciones entre los valores reales y predichos (Fig. S9b–f, paneles de la derecha). Observamos una reducción de la capacidad de generalización cuando se excluyeron los pacientes de mayor edad del entrenamiento (Fig. S9d). Hipotetizamos que esta dificultad en la generalización puede deberse a las propiedades distintas de esta subpoblación de pacientes, que muestra una reducción de la diversidad clonal de las células T. No obstante, los resultados de todos los análisis de subgrupos fueron consistentes, lo que indica que TRIM es robusto a las variaciones demográficas y del tamaño de la muestra.

También evaluamos nuestro modelo en el conjunto de datos de cánceres múltiples ^[11] (Fig. S10a). Para evaluar cuantitativamente las predicciones de TRIM, agrupamos las células en el espacio de ARN utilizando k medias (k = 20) y luego calculamos la diversidad clonal real y predicha y el recuento medio de clones para ese grupo (Fig. S10b). Encontramos un alto grado de correlación en ambos casos, r = 0,97 y r = 0,97, respectivamente. Además, las proyecciones UMAP de los datos reales y los datos predichos de nuestro modelo mostraron que el modelo había aprendido a generar, para los pacientes excluidos, perfiles de expresión de ARN precisos basados en la fuente (tejido normal/tumor) y el tipo de cáncer (cáncer de esófago o ESCA, cáncer de tiroides o THCA, carcinoma endometrial del cuerpo uterino o UCEC) (Fig. S10c, gráficos de la izquierda). Además, el modelo también aprendió a igualar el recuento de pseudoclones en cada región del espacio de ARN con el recuento de clones real correspondiente para esa región (Fig. S10c, gráficos de la derecha). En conjunto, estos resultados indican que el modelo logra una alta precisión tanto a nivel del conjunto de datos como a nivel del grupo.

TRIM predice con precisión la expansión clonal y el estado de activación de las células T después del tratamiento

A continuación, consideramos la resolución más granular: la predicción de la expansión a nivel de células individuales. En el conjunto de datos de CHCC, las muestras de sangre pretratamiento incluyeron 15.147 clones de TCR únicos. De estos, 725 persistieron en la sangre postratamiento, y 185 se expandieron después del tratamiento. Cabe destacar que 66 de estos 185 clones expandidos (36%) también se detectaron en el tumor postratamiento. Y 64 de estos 66 TCR se expandieron en el tumor postratamiento. Esta observación respaldó los hallazgos previos de una conexión entre la dinámica de las células T periféricas (sangre) y las asociadas al tumor ^[19], lo que nos motivó a modelar la expansión de las células T utilizando solo datos derivados de la sangre. Debido al desafío de predecir directamente las secuencias de TCR, como hemos demostrado anteriormente, que las secuencias de ARN de una sola célula y los TCR están relacionadas solo a través del recuento de clones, modelamos la expansión a nivel de clon y restringimos el análisis a los TCR presentes en ambas condiciones. Dado un TCR específico y su recuento de clones en la sangre pretratamiento, nuestro objetivo fue predecir si su recuento de clones aumentaría o disminuiría después del tratamiento (normalizado por el número total de células en cada condición). Para ello, mapeamos el TCR al espacio latente y obtuvimos el vector latente correspondiente z. Luego, reemplazamos la incrustación condicional aprendida para "pretratamiento" y, en su lugar, utilizamos la incrustación aprendida para "postratamiento", y decodificamos como de costumbre. Comparamos el recuento de pseudoclones resultante con el recuento de pseudoclones obtenido decodificando los vectores latentes originales, y si el primero era mayor, predijimos la expansión. Definimos la expansión en los datos reales como un TCR que se observa en ambos puntos temporales y cuyo recuento en la distribución postratamiento es mayor que en la distribución pretratamiento, normalizado por el número total de células en cada condición (consulte Métodos para obtener más detalles). Como se muestra en la Fig. 5c, nuestro modelo logró un área bajo la curva de características operativas del receptor (AUC) de 0,92 para esta tarea de predicción. En contraste, los métodos de referencia, como un clasificador de vecinos más cercanos (KNN) o un perceptrón multicapa (MLP), tuvieron un rendimiento similar al de la línea de base aleatoria (línea diagonal) al predecir una etiqueta binaria de expansión, cuando se les proporcionó como entrada la expresión de ARN, las incrustaciones del autoencoder de TCR o la concatenación de ambos.

De manera similar, evaluamos TRIM en el conjunto de datos de CRC. En la evaluación de exclusión de un solo elemento, TRIM logró un AUC de 0,83 para predecir la expansión clonal postratamiento, mientras que los modelos de referencia (KNN y MLP) tuvieron un rendimiento inferior o cercano al aleatorio (incluso después de una optimización exhaustiva de los hiperparámetros), independientemente de si utilizaron la expresión de ARN, las incrustaciones de TCR o su combinación como entrada (Fig. S11a–d, Métodos). Para evaluar aún más la robustez, probamos TRIM en entornos de exclusión de grupos, donde los pacientes de prueba tenían mediciones postratamiento. TRIM mantuvo un rendimiento competitivo, logrando AUC de 0,68, 0,81 y 0,78 para las evaluaciones de exclusión de mujeres, exclusión de jóvenes y exclusión de más células, respectivamente, mientras que los modelos de referencia nuevamente tuvieron un rendimiento cercano al aleatorio (Fig. S12).

Luego, evaluamos la capacidad de TRIM para predecir la expansión clonal a nivel de célula individual en el conjunto de datos de cánceres múltiples (Fig. S10d). Predijimos si un TCR en particular que se encuentra en el tejido normal adyacente tenía un recuento de clones más alto en el tumor. Realizamos esta tarea de predicción de la misma manera generativa condicional que en el conjunto de datos de CHCC, restringiéndonos solo a los TCR que aparecen en ambas condiciones.

Como antes, los clasificadores KNN y las redes neuronales de referencia que recibieron la concatenación de las representaciones de ARN y TCR tuvieron un rendimiento similar al de la línea de base aleatoria, mientras que nuestro modelo superó claramente a estos modelos de referencia. Sin embargo, el AUC es menor que en el conjunto de datos de CHCC (0,57 frente a 0,92), lo que puede explicarse por el hecho de que las poblaciones de células T en la sangre periférica circundante contienen más información sobre las poblaciones de células T en el tumor que las que se encuentran en el tejido normal adyacente. Hipotetizamos que la capacidad del modelo para predecir la expansión de las células T en el tumor a partir del tejido normal adyacente mejor que al azar se debe a la señalización y la comunicación entre estas poblaciones de células T, aunque en menor medida que con las células T circundantes en la sangre.

Para validar la utilidad de nuestro modelo para los análisis biológicos posteriores, comparamos sus predicciones basadas en los datos de sangre pretratamiento con los análisis de los fenotipos de tumores y células T postratamiento informados en la ref. ^[19] en el conjunto de datos de CHCC (Fig. 5d, e, Fig. S13). En ^[19], los autores informaron sobre la expansión de las células T CD8 proliferativas CD38+/HLA-DR+ en la sangre durante el curso del tratamiento y encontraron que la frecuencia de las células T CD8 activadas KLRG1 negativas PD-1+ en la sangre pretratamiento y en el tratamiento estaba fuertemente asociada con la respuesta patológica intratumoral. Primero, calculamos el porcentaje de células T proliferativas CD38+/HLA-DR+ para cada condición basándonos en nuestros datos predichos y encontramos que esta población aumenta en las condiciones postratamiento (Fig. S13). A nivel de paciente y punto temporal, las predicciones de nuestro modelo están bien correlacionadas con la verdad fundamental para esta estadística (r = 0,89). El parámetro de pendiente de una regresión de los valores reales de esta métrica con los valores predichos inferidos a partir de la salida de nuestro modelo también se muestra, lo que da como resultado β1 = 0,82. A continuación, modelamos la expresión de KLRG1 en las células CD8 PD-1+ y encontramos que disminuye más para los pacientes que tienen una respuesta patológica que para aquellos que no la tienen, con una correlación por paciente de 0,82 (Fig. 5d). Por último, aislamos las células CD8+ predichas, calculamos la proporción que son CD38+/HLA-DR+ y encontramos que esto también está fuertemente correlacionado con la respuesta patológica (Fig. 5e). Las predicciones de nuestro modelo están altamente correlacionadas con los datos reales por paciente (r = 0,90). En general, estos análisis muestran que TRIM genera predicciones precisas y significativas de la expansión clonal.

TRIM identifica genes asociados con la expansión clonal

Para analizar más a fondo el aspecto de predicción de expansión de nuestro modelo, investigamos las características que utiliza para realizar esta tarea. Específicamente, realizamos un análisis de expresión génica diferencial al comparar todos los pseudoclones predichos para expandirse por nuestro modelo con aquellos que no se predijo que se expandieran en las células T de la sangre. Luego, clasificamos los genes según su cambio logarítmico en la expresión, sujeto a la restricción de que se expresen en al menos el 10% de las células. Repetimos este proceso para los datos reales, basándonos en la expansión observada y su expresión génica medida en cada condición. Luego, comparamos las dos listas clasificadas considerando la intersección de sus primeros k elementos, para valores crecientes de k. La superposición de un modelo aleatorio estaría en la línea diagonal, obteniendo solo una intersección de genes en promedio. Primero evaluamos TRIM en el conjunto de datos de CHCC de cabeza y cuello. El AUC de nuestro modelo con respecto a esta línea de base es de 0,92 (Fig. 6a). Repetimos este proceso por separado para los subconjuntos de células T CD4+ y CD8+, lo que resultó en AUC de 0,80 y 0,77, respectivamente (Figura suplementaria 14a). También evaluamos este aspecto de nuestro modelo tomando una firma génica de expansión y una firma génica de no expansión que se identificaron previamente en la referencia ^[19]. Cuando restringimos nuestro análisis a estos genes y los clasificamos por el cambio logarítmico en la expresión, los genes de la firma de expansión tenían más probabilidades de estar sobreexpresados en los clones que el modelo predijo que se expandirían (Fig. 6a). La cuantificación utilizando una prueba de rango de Wilcoxon ^[45] para la hipótesis de si el modelo clasifica los genes de la firma de expansión más alto que los genes de la firma de no expansión resultó en un valor p de 0,064. Dado que las firmas se derivaron de tejido tumoral, mientras que nuestro análisis modela la expansión en la sangre, esperamos cierta atenuación; no obstante, el resultado muestra un enriquecimiento direccionalmente consistente y se alinea con la evidencia previa de que los programas de expansión asociados al tumor son detectables en la sangre. Esto demuestra que nuestro modelo predice con precisión no solo qué clones se expandirán, sino también qué genes se expresarán de manera diferente durante esa expansión.

Llevando este análisis un paso más allá, investigamos el modelo para identificar los genes que se predice que se expresarán de manera diferente (DE) desde el pretratamiento hasta el postratamiento por los pseudoclones predichos para expandirse en la sangre (Figura suplementaria 14b). El análisis de enriquecimiento de conjuntos de genes de los genes DE reveló un enriquecimiento para las vías del ciclo celular y la proliferación celular, que se espera que estén involucradas en la expansión clonal (Figura suplementaria 14c). Estos resultados respaldan la precisión de las predicciones de expresión génica de nuestros pseudoclones.

Luego, evaluamos la capacidad de TRIM para identificar genes expresados de manera diferente en el conjunto de datos de CRC. En la evaluación de exclusión de uno, TRIM logró un AUC de 0,85 para todas las células T, 0,82 para las células CD4 y 0,79 para las células CD8 (Figura suplementaria 15a). En diferentes análisis de subgrupos, TRIM mantuvo un fuerte rendimiento, con AUC para todas las células T de 0,83, 0,90 y 0,81 para la exclusión de mujeres, la exclusión de jóvenes y la exclusión de más células, respectivamente; para las células CD4, 0,83, 0,83 y 0,81; y para las células CD8, 0,72, 0,84 y 0,72 (Figura suplementaria 15b–d). Estos resultados validan aún más la capacidad de TRIM para generalizarse a diferentes tipos de cáncer, así como en diversos grupos de pacientes.

Para identificar los genes que el modelo utiliza para predecir la expansión clonal de las células T CD8, realizamos un análisis de saliencia. Esto implicó la creación de un clasificador para la expansión clonal basado en el espacio latente de la red y el examen de los gradientes en relación con cada gen de entrada. A través de este análisis de saliencia, generamos una clasificación de genes basada en su poder predictivo para la expansión, que va desde aquellos que son fuertemente predictivos negativos hasta aquellos que son fuertemente predictivos positivos. Comparamos esta clasificación (gradiente) con los resultados de los análisis DE de rutina para evaluar su relevancia. Específicamente, examinamos tres escenarios DE: (1) grupos expandidos frente a no expandidos en el postratamiento, (2) grupos que se expanden frente a no expandidos en el pretratamiento y (3) un caso simulado sin datos de postratamiento, comparando células con tamaños de clon grandes frente a tamaños de clon pequeños. Para cada método, identificamos genes importantes basados en valores p corregidos (para los análisis DE) o puntajes de clasificación (en nuestro método). Los análisis de correlación de Pearson mostraron que los genes identificados como importantes para que nuestro modelo prediga la expansión clonal se alinean más estrechamente con la comparación de expansión en el pretratamiento frente a no expansión (0,52) que con los otros dos enfoques DE (-0,07 y 0,22) (Fig. 6b). De la Fig. 6c, observamos que los genes de gradiente tienen una superposición limitada con los genes identificados a través de los cambios logarítmicos en la expresión en las comparaciones de postratamiento y tamaño de clon, y un número relativamente mayor de genes compartidos con la lista DE de pretratamiento.

Para evaluar aún más la relevancia biológica de los genes identificados por cada enfoque, realizamos un análisis de enriquecimiento de términos de Ontología Génica (GO) utilizando el paquete gseapy ^[46] (Figura suplementaria 16, Métodos). Descubrimos que los conjuntos de genes identificados como enriquecidos entre los genes basados en el gradiente coincidieron en gran medida con los conjuntos de genes identificados utilizando los análisis DE más estándar (Figura suplementaria 16). Luego, consultamos a ChatGPT para compilar una lista de genes asociados con la expansión de células T y analizamos las listas resultantes en los diferentes métodos (Fig. 6d). Este análisis se inspiró en los hallazgos de la referencia ^[47], que demostró el potencial de los modelos de lenguaje grandes para anotar conjuntos de genes. Nuestro método ("Gradiente"), junto con los métodos supervisados que utilizaron etiquetas de postratamiento ("Expandido frente a no expandido" y "Expansión frente a no expansión"), identificó con éxito genes involucrados en el ciclo celular (AURKA, CENPU, DUT), la función efectora de las células T y la citotoxicidad (PRF1, GZMA, GZMH, GNLY, NKG7, CST7), la activación y la coestimulación de las células T (CTLA4, TNFRSF9), las vías metabólicas y mitocondriales (ATP5PF, VDAC2, DHTKD1, IBA57) y la migración de las células T (CCR7, CCL5, IL6R, FOXO1). En particular, las células T CD38+ CD8 que se expandieron en la sangre durante el tratamiento se ha demostrado que exhiben firmas de activación y proliferación; aquí, la identificación de una firma de migración de células T es nueva y sugiere que estas células pueden migrar al tumor y a los ganglios linfáticos, lo que es consistente con que esta población tenga TCR que se expandieron en los tumores durante la terapia. Es importante destacar que, si bien los análisis de expresión diferencial (DE) supervisados sirven como límites superiores teóricos, porque se basan directamente en las etiquetas de postratamiento, nuestro método pudo recuperar señales biológicas similares sin requerir mediciones de postratamiento en el momento de la inferencia. Aunque TRIM requiere acceso a datos de pretratamiento y postratamiento emparejados durante el entrenamiento, una vez entrenado, se puede aplicar a nuevas muestras de solo pretratamiento para identificar programas genéticos asociados con la expansión de células T. Esta capacidad proporciona tanto la interpretabilidad del modelo, al demostrar que sus predicciones se basan en señales biológicas significativas, como la relevancia clínica, ya que proporciona estas perspectivas predictivas en un entorno clínico común donde las muestras de postratamiento no están disponibles. Esto permite que TRIM proporcione información significativa en entornos donde las mediciones de postratamiento no están disponibles, cerrando la brecha entre los análisis supervisados y la aplicabilidad en el mundo real. En contraste, los genes expresados de manera diferente identificados al comparar clones grandes frente a clones pequeños en el pretratamiento se enriquecieron en programas citotóxicos, efectores e inmunorreguladores, pero no mostraron un enriquecimiento significativo para los programas relacionados con la expansión (evaluado mediante pruebas hipergeométricas por conjunto de genes con corrección FDR; el enriquecimiento se muestra como valores p ajustados por FDR en la Fig. 6d). Esto sugiere que TRIM proporciona una visión más informativa y biológicamente fundamentada de la expansión de las células T en condiciones del mundo real.

Discusión

Los avances recientes en las tecnologías para perfilar conjuntamente los transcriptomas y las secuencias de TCR a nivel de una sola célula han llevado al desarrollo de métodos computacionales que tienen como objetivo integrar estas modalidades ^[32]–^[34]. Sin embargo, muchos de estos enfoques realizan simplificaciones: al agregar datos (por ejemplo, promediar la expresión génica en las células dentro del mismo clonotipo y descartar la variabilidad a nivel de una sola célula), al restringir sus resultados a predicciones categóricas como etiquetas de clúster o puntajes de especificidad de antígeno, o al proyectar nuevas células en una incrustación aprendida durante el entrenamiento, limitando su capacidad para generar resultados transcriptómicos y clonales completos y generalizarse a entornos fuera de la muestra ^[32]–^[34]. Nuestro enfoque aborda estas limitaciones al modelar la estructura conjunta de la expresión génica y la información de la secuencia de TCR a nivel de una sola célula, sin reducir los datos a promedios por clon. Al modelar explícitamente la identidad clonal y la dinámica de expansión, TRIM permite predicciones que van más allá de la clasificación de etiquetas, generando perfiles de expresión génica completos y comportamientos clonales para células no vistas en momentos posteriores o en diferentes sitios de tejido, superando las limitaciones de los métodos anteriores y desbloqueando nuevas aplicaciones en la inferencia inmunológica.

Utilizamos TRIM para predecir los perfiles transcriptómicos y los clonotipos de las células T en pacientes con CHCC de cabeza y cuello no vistos. Observamos poca relación entre la similitud de la secuencia de TCR de las células y sus valores de expresión génica, en contraste con los informes de que la secuencia de TCR puede predecir los resultados del destino de las células T. Las células T con secuencias de TCR similares pueden compartir programas de transcripción porque sus TCR tienen una especificidad de antígeno, una intensidad de señalización y se encuentran con microambientes de tejido similares, lo que luego impulsa la convergencia hacia estados funcionales similares. Sin embargo, incluso las células T con secuencias de TCR compartidas pueden divergir en sus estados funcionales debido a factores extrínsecos al TCR: exposición variable a antígenos y citocinas, moléculas de coestimulación variables, plasticidad epigenética y bifurcaciones estocásticas del destino celular ^[48],^[49]. Estudios anteriores, incluido CoNGA, metaCoNGA y una variedad de marcos de puntuación de TCR (TiRP, TCR-innate, TCR-CD8, TCR-reg, TCR-mem), examinaron sistemáticamente las relaciones entre la secuencia de TCR y la expresión génica, y establecieron que la secuencia de TCR proporciona información sobre la especificidad del antígeno y el destino de las células ^[32],^[40],^[50],^[51]. Sin embargo, este efecto fue mayor en las células T no convencionales, que no están presentes en nuestros datos (como se muestra por la falta de expresión de los marcadores de células T no convencionales TCRD, TCRG y baja expresión de RORC; Figura suplementaria 17), y los efectos extrínsecos al TCR fueron el determinante más importante del destino de las células T.

A pesar de la débil relación entre la similitud de la secuencia de TCR y la expresión génica, existió una fuerte relación entre la expresión génica y la expansión clonal, lo que permitió una predicción precisa del estado de las células T y la expansión clonal a partir de muestras de sangre de pretratamiento. Mostramos que esta capacidad predictiva se generalizó a otros entornos de enfermedades al evaluar TRIM en el conjunto de datos de CRC y probar la solidez de TRIM en subgrupos de pacientes con diferentes covariables. También validamos que TRIM predice con precisión la expansión clonal a partir de tejidos adyacentes normales en múltiples tipos de cáncer.

Para interpretar mejor nuestro modelo, analizamos qué genes utiliza para predecir la expansión clonal en el estudio de inmunoterapia del CHCC. Nuestro análisis reveló que el modelo se basa en genes con significado biológico, incluidos muchos que se corresponden con marcadores conocidos de la actividad de las células T. En particular, nuestros hallazgos respaldan la idea de que las células T pretratadas, destinadas a expandirse durante el tratamiento, ya presentan un programa transcripcional distinto. Esto coincide con observaciones previas que indican que muchas de las poblaciones de células T en expansión en la sangre durante el bloqueo de puntos de control comparten TCR con las células T infiltrantes en el tumor, lo que sugiere que son reactivas al tumor. Una posible interpretación es que el reconocimiento previo de los antígenos tumorales (posiblemente en los ganglios linfáticos que drenan el tumor) desencadena cambios transcripcionales tempranos relacionados con la activación de las células T, las funciones efectoras y la adaptación metabólica. Esta hipótesis está aún más respaldada por datos que muestran que la frecuencia pretratamiento de las células T CD8+ activadas (CD38+ HLA-DR+) en la sangre se correlaciona fuertemente con la posterior respuesta tumoral al bloqueo de puntos de control.

Más allá de las firmas relacionadas con el sistema inmunitario, nuestro análisis de saliencia también reveló genes involucrados en el metabolismo y las vías de señalización, una categoría que fue menos prominente en los análisis diferenciales de expresión tradicionales. Esto sugiere que la reprogramación metabólica puede ser un determinante temprano de la expansión de las células T, lo que podría influir en la posterior activación inmunitaria. La identificación de genes metabólicos apoya aún más la idea de que las células T en expansión están preparadas antes del tratamiento, posiblemente experimentando cambios metabólicos que permiten la posterior expansión clonal. En conjunto, estos hallazgos resaltan la utilidad de nuestro modelo para descubrir señales biológicamente significativas que informan sobre los mecanismos de activación de las células T, la función efectora y la respuesta terapéutica.

Hemos demostrado que la estructura de nuestro modelo se generaliza bien en diferentes tipos de cáncer^[19],^[27] y es robusta a los factores que introducen variabilidad intra-muestral, como la edad y el sexo, o el número de células secuenciadas por paciente. Si bien reconocemos que un modelo entrenado en múltiples cohortes, que incluya una gran variedad de tipos de cáncer, podría aplicarse directamente a entornos de tipos de cáncer no vistos sin volver a entrenarlo, dado el actual escaso número de conjuntos de datos de este tipo, creemos que esto no es factible en la actualidad, ni es la estrategia más precisa y rentable para abordar el problema que tenemos entre manos: ser capaz de predecir los datos faltantes en cohortes complejas a lo largo del tiempo y en diferentes sitios tisulares. Hasta donde sabemos, ningún otro método computacional utiliza de manera productiva el ARN y las secuencias de TCR para predecir los estados transcripcionales en diferentes tejidos y puntos temporales; sin embargo, el trabajo presentado aquí aún deja muchas preguntas sin respuesta. Por ejemplo, el estado transcripcional de una célula T después de la inmunoterapia depende, en cierta medida, de la unión del antígeno tumoral y, por lo tanto, de la secuencia de TCR. Una pregunta abierta es cómo aprovechar los modelos externos de predicción de la unión de antígenos en este marco.

Existe una necesidad urgente de biomarcadores de la progresión del cáncer y la respuesta al tratamiento, y los modelos multimodales que incorporen mediciones longitudinales en múltiples sitios tisulares impulsarán el descubrimiento de estos biomarcadores. TRIM es el primer modelo en lograr esto, y en general, consideramos que TRIM es una herramienta valiosa para predecir los estados de las células T y la respuesta a las inmunoterapias contra el cáncer en pacientes con cáncer, basándose en mediciones de la sangre periférica. También demostramos que nuestro enfoque se generaliza a las células T en diversos conjuntos de datos de cáncer. Si bien los conjuntos de datos multimodales longitudinales son actualmente limitados, esperamos que se produzcan mejoras en el modelo a medida que se recopilen más conjuntos de datos, lo que impulsará el descubrimiento de programas inmunitarios compartidos que caractericen las respuestas productivas de las células T durante la inmunoterapia y en otros contextos de enfermedades^[19],^[27],^[28].

Métodos

TRIM

Supongamos que los datos se capturan mediante dos matrices, XRNA y XTCR, donde n es el número de células y dRNA y dTCR son las dimensiones de la modalidad de ARN y la modalidad de TCR, respectivamente. Cada célula también tiene etiquetas de condición asociadas, l, que comprenden: paciente lp ∈ {0…npatient}, tiempo de tratamiento lt ∈ {pretratamiento, postratamiento}, origen lb ∈ {sangre, tumor}, donde npatient denota el número de pacientes en el conjunto de datos.

Nuestro modelo consta de dos codificadores, ERNA(x∣l) y ETCR(x∣l), y dos decodificadores, DRNA(z∣l) y DTCR(z∣l). Los codificadores toman como entrada las células de sus respectivas modalidades, xRNA ~ XRNA y xTCR ~ XTCR, y producen representaciones latentes μRNA, σRNA, μTCR, σTCR en el espacio latente, donde dz es la dimensionalidad elegida para el espacio latente. Estas representaciones latentes se promedian en un solo vector a lo largo de las dimensiones latentes, μ = mean(μRNA, μTCR) y σ = mean(σRNA, σTCR). Luego, se utilizan para parametrizar una distribución normal de la cual se muestrea y que luego se proporciona al decodificador junto con las condiciones originales para producir las salidas xRNA y xTCR, respectivamente. En conjunto, esto produce:

Los codificadores y decodificadores se entrenan conjuntamente bajo un régimen de auto-codificador variacional^[52]. Existe una pérdida de divergencia de Kullback-Leibler que anima a que el vector latente se distribuya como una normal:

Pérdida de clonalidad contrastiva

La pérdida de reconstrucción anima a que la salida se acerque a los puntos originales. Para el espacio de ARN, utilizamos la pérdida L2 estándar. Para el espacio de TCR, introducimos la pérdida de clonalidad contrastiva:

donde Iij = 1 si i = j y Iij = 0 si i ≠ j, y δ es un margen elegido. Por lo tanto, la pérdida de reconstrucción completa es:

La pérdida de clonalidad contrastiva opera sobre la salida del decodificador en el espacio de TCR; utiliza una variable indicadora para seleccionar las reconstrucciones de las células que son clones de TCR idénticos y minimiza la distancia entre sus reconstrucciones. Para todos los demás pares de células, sus reconstrucciones en el espacio de TCR deben estar separadas por al menos el margen δ. Por lo tanto, el espacio de TCR reconstruido conserva algunos aspectos del espacio original, a saber, que los clones de TCR están cerca. El espacio de TCR reconstruido no se ve obligado a preservar la estructura global fuera de los clones. La idea de preservar solo la información de los clones locales en nuestro modelo conjunto de ARN y TCR, en lugar de a) preservar la información del espacio de TCR global o b) no preservar ninguna información del espacio de TCR, se basa en nuestras observaciones en las Figuras 3 y 4. Los clones de TCR están más cerca en el espacio de ARN que los aleatorios, incluso después de tener en cuenta factores como el tipo de célula, por lo que incluir esa información es valiosa para la representación. Pero fuera de los clones exactos, ningún incrustamiento del espacio de TCR produjo una estructura de distancia que fuera predictiva del espacio de ARN. En cambio, utilizando la pérdida de reconstrucción L2 estándar, el modelo tendría que almacenar suficiente información en el espacio latente z para reconstruir completamente las representaciones de TCR de todos los puntos. Dado que el espacio de TCR no comparte la misma estructura de distancia por pares que el espacio de ARN, obligar al modelo a reconstruir completamente el espacio de TCR empeoraría su rendimiento en el espacio de ARN. En resumen, nuestra pérdida de clonalidad contrastiva obliga al modelo a aprender una reconstrucción que preserve solo la estructura local de TCR informativa para el espacio de ARN.

Asignación de pseudoclones

Dado que nuestro modelo generativo muestrea a partir de una distribución aleatoria continua, su salida decodificada en el espacio de TCR también es continua. Para comparar métricas basadas en el recuento de clones, que, como hemos demostrado en las Figuras 3 y 4, es la única señal compartida entre los espacios de ARN y TCR, la salida generada debe discretizarse. El procedimiento que utilizamos para hacerlo es un algoritmo voraz que elige la representación de TCR de una célula generada aleatoriamente y asigna todas las demás representaciones de TCR de las células al mismo clon si están dentro de δ (el margen utilizado en la pérdida de clonalidad contrastiva) y pertenecen al mismo paciente. Llamamos a estos clones asignados pseudoclones, enfatizando que el espacio de TCR generado no es globalmente similar al espacio de TCR real. Este proceso continúa hasta que todas las células hayan recibido una asignación de pseudoclones. Durante el preprocesamiento de los datos, los recuentos de clones de TCR se calcularon a nivel de conjunto de datos en lugar de por paciente. Esto puede afectar ligeramente los recuentos de expansión para los TCR que se comparten entre pacientes, aunque estos casos son raros (menos del 1 %).

Aprendizaje de incrustaciones de condiciones

TRIM es un modelo condicional que toma como entrada las condiciones correspondientes al paciente, la fuente y el momento de la célula que se codifica o decodifica. Estas condiciones se representan como vectores de longitud fija que se aprenden conjuntamente con el resto del modelo mediante el uso de una red de incrustación separada, F. Como se muestra en la Figura 2a, esta red toma como entrada las células de un subconjunto fijo de nemb células de entrenamiento, una colección a la que nos referiremos como Semb, y las proyecta en un espacio de dimensión demb. Estas mismas células se agrupan y se promedian para calcular el vector de incrustación para cada etiqueta condicional diferente. Por ejemplo, el vector de dimensión demb para el momento de tiempo de pretratamiento se calcula promediando la representación en el espacio de salida de F para los puntos en Semb que son de pretratamiento, cada uno de los cuales es un vector de dimensión demb. En general, la representación para la condición c se calcula como:

donde xi y li son el vector de secuenciación de ARN de una sola célula y la condición de la célula i, respectivamente. Para la incrustación del paciente, primero se subconjunta el conjunto de puntos de incrustación para incluir solo la fuente de sangre y el momento de tiempo de pretratamiento, ya que nuestro caso de uso requiere la extrapolación a un paciente en el que solo estas células están disponibles. Es importante tener en cuenta que el modelo no utiliza el ID de paciente arbitrario como una característica directa. En cambio, el ID de paciente sirve como una clave de agrupación. El modelo aprende un vector de incrustación específico del paciente agregando los datos de ARN y TCR derivados de la sangre de pretratamiento para ese paciente en particular (como se describe anteriormente). Este vector aprendido, que captura la línea de base biológica única del paciente, es lo que se utiliza como la condición a nivel del paciente.

Estas incrustaciones se aprenden utilizando tanto la representación de ARN como la de TCR de cada célula, es decir, F consta de dos redes, FRNA y FTCR, y F(x) es el promedio por dimensiones de FRNA(xRNA) y FTCR(xTCR).

La motivación detrás de este método de aprendizaje de incrustaciones es que la forma estándar de aprender un vector condicional sería indexar en una matriz de incrustación iniciada aleatoriamente cuyos parámetros se aprenden durante el entrenamiento. Esto funciona para los vectores de fuente y momento de tiempo, porque son variables de baja cardinalidad (cada una puede tomar solo dos valores posibles), y nuestras predicciones no tendrán lugar en valores nuevos de estas variables. Sin embargo, para la incrustación del paciente, no hay una forma obvia y efectiva de aprender un vector apropiado que permita una extrapolación precisa a nuevos pacientes no vistos. Intentamos utilizar los metadatos del paciente, pero el modelo no pudo generar con precisión distribuciones de ARN o TCR informativas. Encontramos que solo la distribución de ARN o TCR del paciente en una fuente/momento de tiempo podría predecir con precisión la distribución de ARN o TCR de ese paciente en otra fuente/momento de tiempo. Por lo tanto, formulamos el problema como una generación condicional de las distribuciones de sangre de pretratamiento, tumor de pretratamiento y tumor de postratamiento de un paciente en función de una incrustación de paciente aprendida a partir de la sangre de pretratamiento de ese paciente y las incrustaciones generales para el efecto de fuente y momento de tiempo aprendidas de esta manera.

Entropía de Shannon

Se utilizó la entropía de Shannon para cuantificar la diversidad de los clones de TCR. Para cada muestra, identificamos todos los clones de células T observados, contamos la frecuencia de cada clon único y convertimos estos recuentos en probabilidades dividiendo por el número total de observaciones. Luego, se calculó la entropía como:

donde pi representa la probabilidad del clon i.

Análisis de saliencia

Se realizó un análisis de relevancia del conjunto de datos de CHCC para evaluar cómo la expansión de los TCR en las células T CD8 se predice a partir de los datos sanguíneos previos al tratamiento dentro del marco del modelo presentado. Para ello, se reemplazó la capa latente del autoencoder con una capa de predicción binaria que genera un vector "one-hot" que corresponde a si el TCR de la célula de entrada se expande en la muestra posterior al tratamiento, sin utilizar la capa decodificadora. Se entrenó una versión del modelo que era idéntica en todos los aspectos y se utilizó la función autograd de PyTorch para evaluar el gradiente de la salida de la predicción de expansión con respecto a cada gen de entrada ^[53]. Los gradientes se promediaron luego en todas las células de la distribución para obtener la media informada.

Análisis de la expresión génica diferencial

Para comparar con el análisis de rutina, se calcularon los cambios logarítmicos en la expresión para los genes bajo las siguientes condiciones: entre los clones expandidos y no expandidos después del tratamiento; entre los clones en expansión y los que no lo están, antes del tratamiento; entre los clones con tamaños mayores de 10 y menores de 10 antes del tratamiento; y entre los clones con tamaños mayores de 1 y menores de 1 antes del tratamiento. Estos cambios logarítmicos se calcularon utilizando la función FindMarkers del paquete Seurat ^[54].

Para identificar los genes significativos detectados por cada método, se aplicaron criterios específicos: para el método del cambio logarítmico, se seleccionaron los genes con un valor p ajustado < 0,01 y, para el método basado en el gradiente, se seleccionó el 1% superior de los genes con los gradientes más altos o más bajos. Para comprender mejor la relevancia biológica de estos genes, se realizó un análisis de enriquecimiento de conjuntos de genes utilizando el paquete GSEApy ^[46]. Para la visualización, se seleccionaron los 15 términos GO (valor p ajustado < 0,1) más significativos enriquecidos de cada método y se representó la unión de estos términos en la Figura Suplementaria 16. Los resultados completos del enriquecimiento se proporcionan en los Datos Suplementarios 1. Además, se compiló la unión de los genes significativos identificados en todos los métodos y se utilizó ChatGPT para identificar y resaltar los conjuntos de genes específicamente relevantes para las expansiones de células T, basándose en los hallazgos de la referencia ^[47].

Conjuntos de datos

Conjunto de datos de CHCC

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[19]. Se utilizaron datos emparejados de secuenciación de ARN de célula única (scRNAseq) y secuenciación de TCR de célula única (scTCRseq) de 28 pacientes en 4 condiciones diferentes: sangre antes del tratamiento, sangre después del tratamiento, tumor antes del tratamiento y tumor después del tratamiento. Los puntos de tiempo de antes y después del tratamiento corresponden a los puntos de tiempo B1 y B2 en el artículo de referencia, respectivamente.

Conjunto de datos de cáncer colorrectal

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[27]. Se utilizaron datos emparejados de scRNAseq y scTCRseq de 22 pacientes en 4 condiciones diferentes: sangre antes del tratamiento, sangre después del tratamiento, tumor antes del tratamiento y tumor después del tratamiento. En el estudio original, "antes" y "después del tratamiento" corresponden a las visitas de referencia y durante el tratamiento, respectivamente. Dado que algunos pacientes tenían múltiples muestras durante el tratamiento, se agregaron esas mediciones a nivel del paciente para cada condición para aprovechar al máximo los datos.

Conjunto de datos pan-cáncer

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[11]. Se aplicó un filtro para los pacientes con muestras que contenían al menos 1000 células tanto en el tejido tumoral como en el tejido normal adyacente. Esto resultó en 19 pacientes y 3 tipos de cáncer: ESCA (cáncer de esófago, 5 pacientes), THCA (carcinoma de tiroides, 8 pacientes) y UCEC (carcinoma endometrial del cuerpo uterino, 6 pacientes).

Procesamiento de datos de célula única

Se transformaron los datos de expresión de ARN de cada célula en escala logarítmica y se normalizaron, y se entrenó nuestro modelo en los primeros 100 componentes principales de los datos. Se rellenó la secuencia de TCR hasta una longitud de 80 (la longitud máxima) de las representaciones "one-hot" de los aminoácidos y se añadió un token nulo al realizar el procesamiento de cadenas. Se normalizaron los recuentos de TCR dividiéndolos por el número total de células en cada porción de interés (por ejemplo, al comparar los recuentos previos al tratamiento con los recuentos posteriores al tratamiento, se utilizó el número total de células correspondiente en cada distribución como denominador). Se hizo lo mismo para convertir los recuentos de pseudoclones generados en proporciones para la comparación. La expansión en una variable, por ejemplo, "tratamiento", se definió entonces como un TCR que se observa en ambas distribuciones y que tiene una proporción mayor en las distribuciones posteriores al tratamiento que en las distribuciones previas al tratamiento.

Representaciones de TCR

Las representaciones de TCR utilizadas en TRIM se basan en las incrustaciones entrenadas producidas por un autoencoder CNN entrenado en un paso de preprocesamiento. Cada TCR se representó mediante su secuencia CDR3 de la cadena β, que se introdujo en el modelo como una secuencia de aminoácidos rellenada con nulos. Cada aminoácido se codificó utilizando sus factores de Atchley correspondientes. El autoencoder CNN consta de un codificador con 3 capas, normalización por lotes y activaciones ReLU con fugas que se mapean a un espacio latente de 100 dimensiones. El decodificador es una red convolucional simétrica. Se entrena durante 50.000 pasos y, a continuación, se extrae la capa latente intermedia y se utiliza como entrada para TRIM.

Predicción del paciente excluido

Para maximizar la utilidad de los datos, se entrenó con los datos de todos los pacientes, excepto uno, y se utilizó el paciente excluido para las pruebas. Se repitió esto en los pacientes en un proceso de validación cruzada para obtener el conjunto completo de predicciones. Para evaluar la robustez del modelo, también se realizó una evaluación de exclusión de grupos, en la que se excluyeron los pacientes con características predefinidas del entrenamiento y se evaluó el modelo en todo el grupo excluido.

Al generar células para cada paciente, se utilizó un proceso de muestreo totalmente generativo que no se limita a auto-codificar/decodificar los puntos de entrada reales. En concreto, se muestreó aleatoriamente a partir de una gaussiana isotrópica, se concatenó el punto de tiempo objetivo y la incrustación de origen (por ejemplo, "tumor" y "post-tratamiento") y, a continuación, se decodificó: Como resultado, las predicciones de nuestro modelo generativo no se pueden comparar con los datos reales a nivel de célula, ya que cada z muestreado de la distribución gaussiana es aleatorio. Sin embargo, la distribución de células por condición se puede comparar con la distribución real de células para esa condición, y lo mismo ocurre con cualquier otro subconjunto identificable de condiciones, o incluso con la distribución en todo el conjunto de datos.

Métrica para la evaluación comparativa

El rendimiento de TRIM y los métodos de referencia para predecir la expansión de TCR se evaluó utilizando la métrica del área bajo la curva de características operativas del receptor (ROC) (AUC). La AUC mide la capacidad del modelo para distinguir entre clases positivas y negativas calculando el área bajo la curva ROC, que representa la tasa de verdaderos positivos frente a la tasa de falsos positivos en diferentes umbrales de clasificación. Esta métrica evalúa la capacidad discriminatoria general de un modelo; los valores de AUC más altos indican un mejor rendimiento.

Resumen del informe

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

TRIM

Sean los datos capturados por dos matrices y donde n es el número de células y dRNA y dTCR son las dimensiones de la modalidad de ARN y la modalidad de TCR, respectivamente. Cada célula también tiene etiquetas de condición asociadas l compuestas por: paciente lp ∈ {0…npaciente}, tiempo de tratamiento lt ∈ {previo al tratamiento, posterior al tratamiento}, origen lb ∈ {sangre, tumor}, donde npaciente denota el número de pacientes en el conjunto de datos.

Nuestro modelo consta de dos codificadores ERNA(x∣l), ETCR(x∣l) y dos decodificadores DRNA(z∣l), DTCR(z∣l). Los codificadores toman como entrada células de sus respectivas modalidades xRNA ~ XRNA y xTCR ~ XTCR y producen representaciones latentes μRNA, σRNA, μTCR, σTCR en el espacio latente , donde dz es la dimensionalidad elegida para el espacio latente. Estas representaciones latentes se promedian en un vector a lo largo de las dimensiones latentes μ = media(μRNA, μTCR) y σ = media(σRNA, σTCR). A continuación, se utilizan para parametrizar una distribución normal de la que se muestrea y que luego se proporciona al decodificador junto con las condiciones originales para producir salidas respectivamente. En conjunto, esto produce:

Los codificadores y decodificadores se entrenan conjuntamente en un régimen de autoencoder variacional ^[52]. Existe una pérdida de divergencia KL que anima a que el vector latente se distribuya como una normal:

Pérdida de clonalidad contrastiva

La pérdida de reconstrucción anima a que la salida se acerque a los puntos originales. Para el espacio de ARN, utilizamos la pérdida L2 estándar. Para el espacio de TCR, introducimos la pérdida de clonalidad contrastiva: donde , Iij = 1 si i = j y Iij = 0 si i ≠ j, y δ es un margen elegido. Por lo tanto, la pérdida de reconstrucción completa es:

La pérdida de clonalidad contrastiva opera sobre la salida del decodificador en el espacio de TCR; utiliza una variable indicadora para seleccionar las reconstrucciones de las células que son clones de TCR idénticos y minimiza la distancia entre sus reconstrucciones. Para todos los demás pares de células, sus reconstrucciones en el espacio de TCR deben estar separadas por al menos el margen δ. El espacio de TCR reconstruido conserva así algunos aspectos del espacio original, a saber, que los clones de TCR están cerca. El espacio de TCR reconstruido no se ve obligado a preservar la estructura global fuera de los clones. La idea de preservar solo la información de los clones locales en nuestro modelo conjunto de ARN y TCR, en lugar de a) preservar la información del espacio de TCR global, o b) no preservar ninguna información del espacio de TCR, se basa en nuestras observaciones en las Figuras 3 y 4. Los clones de TCR están más cerca en el espacio de ARN que al azar, incluso después de tener en cuenta factores como el tipo de célula, por lo que incluir esa información es valioso para la representación. Pero fuera de los clones exactos, ninguna incrustación del espacio de TCR produjo una estructura de distancia que fuera predictiva del espacio de ARN. Al utilizar en su lugar la pérdida de reconstrucción L2 estándar, el modelo tendría que almacenar suficiente información en el espacio latente z para reconstruir completamente las representaciones de TCR de todos los puntos. Dado que el espacio de TCR no comparte la misma estructura de distancia por pares que el espacio de ARN, obligar al modelo a reconstruir completamente el espacio de TCR empeoraría su rendimiento en el espacio de ARN. En resumen, nuestra pérdida de clonalidad contrastiva obliga al modelo a aprender una reconstrucción que preserve solo la estructura local de TCR informativa para el espacio de ARN.

Asignación de pseudoclones

Dado que nuestro modelo generativo muestrea a partir de una distribución aleatoria continua, su salida decodificada en el espacio de TCR también es continua. Para comparar las métricas basadas en el recuento de clones, que, como hemos demostrado en las Figuras 3 y 4, es la única señal compartida entre los espacios de ARN y TCR, la salida generada debe discretizarse. El procedimiento que utilizamos para hacerlo es un algoritmo voraz que elige una representación de TCR de una célula generada aleatoriamente y asigna todas las demás representaciones de TCR de las células al mismo clon si están dentro de δ (el margen utilizado en la pérdida de clonalidad contrastiva) y pertenecen al mismo paciente. Llamamos a estos clones asignados pseudoclones, enfatizando que el espacio de TCR generado no es globalmente similar al espacio de TCR real. Este proceso continúa hasta que todas las células han recibido una asignación de pseudoclones. Durante el preprocesamiento de los datos, los recuentos de clones de TCR se calcularon a nivel del conjunto de datos en lugar de por paciente. Esto puede afectar ligeramente los recuentos de expansión de los TCR que se comparten entre pacientes, aunque estos casos son raros (menos del 1%).

Aprendizaje de incrustaciones de condición

TRIM es un modelo condicional que toma como entrada las condiciones correspondientes al paciente, la fuente y el punto temporal de la célula que se está codificando o decodificando. Estas condiciones se representan como vectores de longitud fija que se aprenden conjuntamente con el resto del modelo mediante el uso de una red de incrustación separada (F). Como se muestra en la Figura 2a, esta red recibe células de un subconjunto fijo de nemb células de entrenamiento y las transforma en un espacio de demb dimensiones. Estas mismas células se agrupan y se promedian según cada etiqueta condicional diferente para calcular el vector de incrustación para esa etiqueta condicional. Por ejemplo, el vector de demb dimensiones para el punto temporal pretratamiento se calcula promediando las representaciones en el espacio de salida de F para los puntos en Semb que son pretratamiento, cada uno de los cuales es un vector de demb dimensiones. De forma más general, la representación para la condición c se calcula como: donde xi y li son el vector de secuenciación de ARN de célula única (scRNAseq) y la condición de la célula i, respectivamente. Para la incrustación del paciente, primero seleccionamos los puntos de incrustación para incluir solo la fuente de sangre y el punto temporal pretratamiento, ya que nuestro caso de uso requiere la extrapolación a un paciente en el que solo estas células están disponibles. Es importante señalar que el modelo no utiliza el ID de paciente arbitrario como una característica directa. En cambio, el ID de paciente sirve como una clave de agrupación. El modelo aprende un vector de incrustación específico del paciente mediante la agregación de los datos de ARN y TCR derivados de la sangre pretratamiento para ese paciente específico (como se describe anteriormente). Este vector aprendido, que captura la línea de base biológica única del paciente, es lo que se utiliza como la condición a nivel del paciente.

Estas incrustaciones se aprenden utilizando tanto la representación de ARN como la de TCR de cada célula, es decir, F consta de dos redes, FRNA y FTCR, y F(x) es el promedio por dimensiones de FRNA(xRNA) y FTCR(xTCR).

La motivación detrás de este método de aprendizaje de incrustaciones es que la forma estándar de aprender un vector condicional sería indexar en una matriz de incrustación iniciada aleatoriamente cuyos parámetros se aprenden durante el entrenamiento. Esto funciona para los vectores de fuente y punto temporal, porque estas son variables de baja cardinalidad (cada una puede tomar solo dos valores posibles), y nuestras predicciones no se realizarán en valores nuevos de estas variables. Sin embargo, para la incrustación del paciente, no existe una forma obvia y efectiva de aprender un vector apropiado que permita una extrapolación precisa a nuevos pacientes no vistos. Intentamos utilizar metadatos del paciente, pero el modelo no pudo generar con precisión distribuciones de ARN o TCR informativas. Descubrimos que solo la distribución de ARN o TCR del paciente en una fuente/punto temporal podía predecir con precisión la distribución de ARN o TCR de ese paciente en otra fuente/punto temporal. Por lo tanto, formulamos el problema como una generación condicional de las distribuciones de sangre postratamiento, tumor pretratamiento y tumor postratamiento de un paciente en función de una incrustación del paciente aprendida a partir de la sangre pretratamiento de ese paciente y las incrustaciones generales para el efecto de la fuente y el punto temporal aprendidas de esta manera.

Entropía de Shannon

Se utilizó la entropía de Shannon para cuantificar la diversidad de los clones de TCR. Para cada muestra, identificamos todos los clones de células T observados, contamos la frecuencia de cada clon único y convertimos estos recuentos en probabilidades dividiendo por el número total de observaciones. Luego, la entropía se calculó como donde pi representa la probabilidad del clon i.

Análisis de saliencia

Se realizó el análisis de saliencia del conjunto de datos de CHCC para evaluar cómo se predice la expansión de TCR en las células T CD8 a partir de los datos de sangre pretratamiento dentro del marco del modelo presentado. Para hacer esto, reemplazamos la capa latente del autoencoder con una capa de predicción binaria que genera un vector de un solo punto correspondiente a si el TCR de la célula de entrada se expande en la muestra postratamiento, sin utilizar la capa del decodificador. Entrenamos una versión del modelo que es, de otro modo, idéntica y utilizamos la función autograd de PyTorch para evaluar el gradiente de la salida de la predicción de expansión con respecto a cada gen de entrada^[53]. Luego, los gradientes se promediaron en todas las células de la distribución para obtener la media informada.

Análisis de la expresión génica diferencial

Para comparar con el análisis de rutina, calculamos los cambios logarítmicos en la expresión para los genes bajo las siguientes condiciones: entre los clones expandidos y no expandidos postratamiento; entre los clones en expansión y los no en expansión pretratamiento; entre los clones con tamaños mayores de 10 y menores de 10 pretratamiento; y entre los clones con tamaños mayores de 1 y menores de 1 pretratamiento. Estos cambios logarítmicos se calcularon utilizando la función FindMarkers del paquete Seurat^[54].

Para identificar los genes significativos detectados por cada método, aplicamos criterios específicos: para el método del cambio logarítmico, seleccionamos los genes con un valor p ajustado < 0,01 y, para el método basado en el gradiente, seleccionamos el 1% superior de los genes con los gradientes más altos o más bajos. Para comprender mejor la relevancia biológica de estos genes, realizamos un análisis de enriquecimiento de conjuntos de genes utilizando el paquete GSEApy^[46]. Para la visualización, seleccionamos los 15 términos GO significativos más importantes (valor p ajustado < 0,1) de cada método y trazamos la unión de estos términos en la Figura Suplementaria 16. Los resultados completos del enriquecimiento se proporcionan en los Datos Suplementarios 1. Además, compilamos la unión de los genes significativos identificados en todos los métodos y utilizamos ChatGPT para identificar y resaltar los conjuntos de genes específicamente relevantes para la expansión de las células T, basándonos en los hallazgos de la referencia ^[47].

Conjuntos de datos

Conjunto de datos de CHCC

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[19]. Utilizamos datos de scRNAseq y scTCRseq emparejados de 28 pacientes en 4 condiciones diferentes: sangre pretratamiento, sangre postratamiento, tumor pretratamiento y tumor postratamiento. Los puntos temporales pretratamiento y postratamiento son los puntos temporales B1 y B2 en el artículo de referencia, respectivamente.

Conjunto de datos de cáncer colorrectal

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[27]. Utilizamos datos de scRNAseq y scTCRseq emparejados de 22 pacientes en 4 condiciones diferentes: sangre pretratamiento, sangre postratamiento, tumor pretratamiento y tumor postratamiento. En el estudio de origen, "pre" y "postratamiento" corresponden a las visitas de referencia y durante el tratamiento. Debido a que algunos pacientes tenían múltiples muestras durante el tratamiento, agregamos esas mediciones a nivel del paciente para cada condición para aprovechar al máximo los datos.

Conjunto de datos pancanceroso

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[11]. Aplicamos un filtro para los pacientes con muestras que contienen al menos 1000 células tanto en el tumor como en el tejido normal adyacente. Esto resultó en 19 pacientes y 3 tipos de cáncer: ESCA (cáncer de esófago, 5 pacientes), THCA (carcinoma de tiroides, 8 pacientes) y UCEC (carcinoma endometrial del cuerpo uterino, 6 pacientes).

Conjunto de datos de CHCC

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[19]. Utilizamos datos de scRNAseq y scTCRseq emparejados de 28 pacientes en 4 condiciones diferentes: sangre pretratamiento, sangre postratamiento, tumor pretratamiento y tumor postratamiento. Los puntos temporales pretratamiento y postratamiento son los puntos temporales B1 y B2 en el artículo de referencia, respectivamente.

Conjunto de datos de cáncer colorrectal

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[27]. Utilizamos datos de scRNAseq y scTCRseq emparejados de 22 pacientes en 4 condiciones diferentes: sangre pretratamiento, sangre postratamiento, tumor pretratamiento y tumor postratamiento. En el estudio de origen, "pre" y "postratamiento" corresponden a las visitas de referencia y durante el tratamiento. Debido a que algunos pacientes tenían múltiples muestras durante el tratamiento, agregamos esas mediciones a nivel del paciente para cada condición para aprovechar al máximo los datos.

Conjunto de datos pancanceroso

Este conjunto de datos se deriva de los datos disponibles públicamente de la referencia ^[11]. Aplicamos un filtro para los pacientes con muestras que contienen al menos 1000 células tanto en el tumor como en el tejido normal adyacente. Esto resultó en 19 pacientes y 3 tipos de cáncer: ESCA (cáncer de esófago, 5 pacientes), THCA (carcinoma de tiroides, 8 pacientes) y UCEC (carcinoma endometrial del cuerpo uterino, 6 pacientes).

Procesamiento de datos de células únicas

Transformamos los datos de expresión de ARN de cada célula en escala logarítmica y los normalizamos, y entrenamos nuestro modelo en los primeros 100 componentes principales de los datos. Rellenamos la secuencia de TCR a una longitud de 80 (la longitud máxima) de las representaciones de un solo punto de los aminoácidos y un token nulo al realizar el procesamiento de cadenas. Normalizamos los recuentos de TCR por el número total de células en cada porción de interés (por ejemplo, al comparar los recuentos pretratamiento con los recuentos postratamiento, el número total de células correspondiente en cada distribución se utilizó como denominador). Lo mismo se hizo para convertir los recuentos de pseudoclones generados en proporciones para la comparación. La expansión en una variable, por ejemplo, "tratamiento", se definió entonces como un TCR que se observa en ambas distribuciones y que tiene una proporción más alta en la distribución postratamiento que en la distribución pretratamiento.

Representaciones de TCR

Las representaciones de TCR utilizadas en TRIM se basan en las incrustaciones entrenadas producidas por un autoencoder CNN entrenado en un paso de preprocesamiento. Cada TCR se representó mediante su secuencia CDR3 de la cadena β, que se introdujo en el modelo como una secuencia de aminoácidos rellena con nulos. Cada aminoácido se codificó utilizando sus factores de Atchley correspondientes. El autoencoder CNN consta de un codificador con 3 capas, normalización por lotes y activaciones ReLU con fugas que se mapean a un espacio latente de 100 dimensiones. El decodificador es una red convolucional simétrica. Se entrena durante 50.000 pasos y, a continuación, se extrae la capa latente central y se utiliza como entrada para TRIM.

Predicción del paciente excluido

Para maximizar la utilidad de los datos, entrenamos con todos los datos excepto los de un paciente y utilizamos el paciente excluido para las pruebas. Repetimos esto en los pacientes en un proceso de validación cruzada para obtener el conjunto completo de predicciones. Para evaluar la solidez del modelo, también realizamos una evaluación de exclusión de grupos, en la que los pacientes con características predefinidas se excluyeron del entrenamiento y el modelo se evaluó en todo el grupo excluido.

Al generar células para cada paciente, utilizamos un proceso de muestreo totalmente generativo que no simplemente vuelve a codificar/decodificar los puntos de entrada reales. En concreto, muestreamos aleatoriamente a partir de una gaussiana isotrópica, concatenamos la incrustación del punto temporal y la fuente objetivo (por ejemplo, "tumor" y "postratamiento") y, a continuación, decodificamos: Como resultado, las predicciones de nuestro modelo generativo no se pueden comparar con los datos reales a nivel de célula, ya que cada z muestreado de la distribución gaussiana es aleatorio. Sin embargo, la distribución de células por condición se puede comparar con la distribución real de células para esa condición, y lo mismo ocurre con cualquier otro subconjunto identificable de condiciones, o incluso con la distribución en todo el conjunto de datos.

Métrica para la evaluación comparativa

El rendimiento de TRIM y los métodos de referencia para predecir la expansión de TCR se evaluó utilizando la métrica del área bajo la curva de características operativas del receptor (ROC) (AUC). La AUC mide la capacidad del modelo para distinguir entre las clases positivas y negativas calculando el área bajo la curva ROC, que traza la tasa de verdaderos positivos frente a la tasa de falsos positivos en diferentes umbrales de clasificación. Esta métrica evalúa la capacidad discriminatoria general de un modelo; los valores de AUC más altos indican un mejor rendimiento.

Resumen de la información

Se proporciona 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.