L'étude a été menée conformément à la Déclaration d'Helsinki, et le protocole a été approuvé par le Comité d'éthique de l'Hôpital thoracique d'Anhui (K2025-007) le 22 avril 2025. Le consentement éclairé a été obtenu de tous les sujets ayant participé à l'étude.
Extraction et normalisation des données
Les profils transcriptomiques et les jeux de données cliniques correspondants pour le CUPA ont été obtenus à partir des cohortes TCGA et GEO. Le jeu de données TCGA-CUPA a été désigné comme ensemble d'apprentissage, tandis que GSE72094, GSE31210 et GSE26939 ont servi de cohortes pour la validation externe (Tableau 1). De plus, 900 gènes liés au MCR ont été recueillis à partir d'une étude précédente12 (Tableau supplémentaire 1). Les données transcriptomiques ont été annotées à l'aide de GENCODE v36 ou des fichiers d'annotation GPL correspondants. Les identifiants de sondes ont été convertis en symboles géniques, les gènes en double ont été fusionnés à l'aide de la fonction avereps, et seuls les gènes codant pour des protéines ont été conservés afin de générer des matrices d'expression au niveau génique. Pour l'ensemble d'apprentissage TCGA-CUPA, les gènes dont les Fragments par kilobase du modèle d'exon par million de fragments cartographiés (FPKM) étaient < 1 dans plus de 50 % des échantillons ont été éliminés, et les valeurs d'expression restantes ont été transformées en log2 (log2[FPKM+1]). Pour les cohortes de validation GEO, les données d'expression brutes ont été téléchargées, les identifiants de sondes ont été associés aux symboles géniques à l'aide des fichiers d'annotation propres à chaque plateforme, et les sondes multiples correspondant au même gène ont été regroupées en faisant la moyenne de leurs valeurs d'expression. Ces jeux de données ont été transformés en log2 si nécessaire. Aucune correction d'effet de lot inter-plateformes n'a été appliquée entre TCGA et GEO, car nous avons adopté une stratégie de standardisation par cohorte afin d'assurer une comparabilité relative. Plus précisément, pour les cohortes d'apprentissage et de validation, les valeurs d'expression génique ont été centrées et normalisées (transformation en score z) à l'aide de la moyenne et de l'écart-type de chaque jeu de données individuellement. Les mêmes coefficients de régression de Cox dérivés de l'ensemble d'apprentissage ont ensuite été utilisés pour calculer les scores de risque pour toutes les cohortes. Afin de préserver l'applicabilité clinique et d'éviter un surajustement à l'une quelconque des cohortes de validation, la médiane du score de risque de la cohorte d'apprentissage a été utilisée comme seuil fixe pour répartir les patients en groupes à haut et à faible risque à travers toutes les cohortes de validation externe. Les informations cliniques, incluant l'âge, le sexe, le stade pathologique, le stade Tumeur-Node-Métastase (TNM), le type histologique, la durée de survie, le statut de survie et le type de tissu, ont été extraites lorsque disponibles. Le critère de jugement était la survie globale (SG). Les échantillons présentant des informations de survie incomplètes ou une durée de survie < 30 jours ont été exclus. La durée de survie a été convertie en années, et le statut de survie a été codé 0 pour vivant et 1 pour décédé.
Identification et analyses fonctionnelles de gènes candidats
Le package Limma a identifié les gènes différentiellement exprimés (DEGs) entre les échantillons tumoraux et normaux de LUAD dans le jeu d'apprentissage13. Les critères suivants ont défini les DEGs : |log2FC| > 0,5 et valeur p ajustée < 0,05. Par la suite, l'algorithme de clustering flou mfuzz du package R ClusterGVis a été utilisé pour répartir les DEGs en clusters d'expression distincts. Une analyse de l'ontologie génique – processus biologique (GO-BP) a été réalisée sur les cinq gènes les plus représentatifs de chaque cluster, en fonction de leurs scores d'appartenance. Un ensemble de gènes communs a été obtenu en prenant l'intersection des DEGs avec les MCRGs. Une analyse d'enrichissement fonctionnel utilisant l'ontologie génique/encyclopédie de Kyoto des gènes et des génomes (GO/KEGG) a évalué la pertinence biologique des gènes en chevauchement. Les réseaux d'interactions protéine-protéine (PPI) provenaient de la base de données STRING14. Seules les interactions ayant des scores de confiance > 0,7 ont été conservées afin d'améliorer la fiabilité du réseau.
Recherche de gènes pronostiques
Le package Survival a été utilisé pour réaliser une analyse de régression de Cox univariée afin d'identifier les gènes probablement associés au taux de survie global dans le LUAD15. Les gènes ayant une valeur p < 0,05 ont été considérés comme des indicateurs pronostiques potentiels. La cohorte d'apprentissage TCGA-LUAD comprenait 500 patients disposant de données complètes sur la survie, parmi lesquels 216 (43,2 %) avaient présenté un décès pendant le suivi. Le rapport entre le nombre de gènes candidats (n = 108) et le nombre d'événements (n = 216) était d'environ 1:2, ce qui est acceptable pour une analyse de régression de Cox. Ensuite, une analyse de régression par sélection et réduction absolue minimale (LASSO) ainsi qu'un modèle Extreme Gradient Boosting (XGBoost) ont permis une sélection supplémentaire de caractéristiques. Des modèles de risques proportionnels de Cox ont été construits avec family = "cox" à l'aide de la fonction cv.glmnet du package glmnet. Le paramètre de régularisation optimal a été déterminé par une validation croisée en 10 parties, la valeur λ.min représentant l'erreur minimale de validation croisée, choisie comme valeur optimale de λ. Les gènes présentant des coefficients de régression non nuls ont été extraits comme caractéristiques candidates. Pour le modèle XGBoost, le temps de survie et le statut de survie ont été combinés en une seule variable réponse, des valeurs positives étant attribuées aux décès et des valeurs négatives aux cas censurés. Les paramètres ont été définis comme objective = "survival: cox" et eval_metric = "cox-nloglik", avec 100 itérations et un taux d'apprentissage de 0,1. Après l'entraînement du modèle, les scores d'importance des gènes ont été calculés à partir des gains des caractéristiques. Les 20 gènes les plus importants ont été conservés après tri des scores d'importance par ordre décroissant, afin de réduire la dimensionnalité des caractéristiques et la complexité du modèle. Les gènes communs aux résultats de LASSO et de XGBoost ont été identifiés comme des gènes pronostiques candidats.
Construction et évaluation d'un modèle pronostique
Un modèle pronostique a été développé à l'aide d'une analyse de régression de Cox multivariée portant sur les gènes candidats identifiés. Les scores de risque ont été calculés individuellement comme suit :
.
où Coefi désigne le coefficient pour le gène i, et Expi indique la valeur d'expression génique respective. Les individus ont ensuite été répartis en deux groupes : à haut risque et à faible risque, en utilisant le score médian de risque comme seuil de division. Ensuite, des courbes caractéristiques de fonctionnement du récepteur (ROC) dépendant du temps ont été établies. Afin d'évaluer le risque de surajustement, une validation interne par rééchantillonnage bootstrap comprenant 1 000 itérations a été réalisée afin de calculer l'indice C corrigé du biais et les AUC dépendant du temps, avec des intervalles de confiance à 95 %. Des courbes de calibration ont été générées pour évaluer la concordance entre les probabilités de survie prédites et observées à 2, 3 et 5 ans. De plus, une analyse de courbe décisionnelle (DCA) a été effectuée à l'aide du package ggDCA dans R afin d'évaluer le bénéfice net clinique du modèle aux points temporels de 2, 3 et 5 ans, quantifiant ainsi la valeur potentielle du score de risque dans la prise de décision clinique selon différentes probabilités seuils. Les différences de survie entre les groupes stratifiés par risque et selon d'autres catégories cliniques ont été comparées à l'aide de courbes de survie de Kaplan-Meier (KM) avec un test du log-rank. En outre, l'analyse Shapley Additive exPlanations (SHAP) a été utilisée afin d'élucider la contribution individuelle des gènes à la performance du modèle, permettant une interprétation explicative a posteriori.
Développement et validation externe d'un nomogramme
Les relations entre les scores de risque calculés et diverses caractéristiques cliniques (y compris le sexe, l'âge et le stade TNM) ont été examinées à l'aide de tests de Wilcoxon ou de Kruskal-Wallis afin d'évaluer l'applicabilité clinique du modèle. Pour déterminer si le score de risque constituait un facteur pronostique indépendant, les variables cliniques ainsi que le score de risque ont été intégrées dans un modèle de régression de Cox multivariée. Ensuite, un nomogramme pronostique combinant les facteurs de risque cliniques indépendants (par exemple, le stade) et le score de risque génétique a été construit à l'aide du package R regplot afin de personnaliser les prédictions de probabilité de survie. Des courbes de calibration ont été utilisées pour évaluer la concordance entre la probabilité de survie prédite par le nomogramme et les résultats de survie réels. Enfin, la capacité prédictive finale et la généralisabilité du système de nomogramme intégré ont été rigoureusement validées à l'aide de courbes ROC dépendantes du temps et d'analyses exhaustives de sous-groupes cliniques KM à travers les cohortes.
Analyses de l'infiltration immunitaire et des sous-types immunitaires
CIBERSORT, en utilisant la matrice de signatures des gènes des leucocytes (LM22), a été utilisé pour estimer les proportions relatives des 22 types de cellules immunitaires afin d'évaluer l'infiltration des cellules immunitaires chez les patients atteints de LUAD. Les relations entre les niveaux d'expression des gènes pronostiques et l'infiltration immunologique ont été évaluées par une analyse de corrélation de Spearman. Les scores immunitaire, stromal, de pureté tumorale et d'estimation ont été obtenus à l'aide de l'algorithme ESTIMATE, et le test de Wilcoxon a été utilisé pour évaluer les différences entre les groupes à risque. Les patients atteints de LUAD ont été répartis en six sous-types immunitaires à l'aide du package ImmuneSubtypeClassifier16. Le test de Wilcoxon a également été utilisé pour comparer les répartitions des sous-types immunitaires entre les groupes à risque.
Analyses des points de contrôle immunitaires, du score d'immunophénotype et du cycle de l'immunité anticancéreuse
Dans cette étude, le test de Wilcoxon à deux échantillons indépendants a été utilisé pour évaluer 21 gènes de points de contrôle immunitaire17 dans des groupes stratifiés selon le risque, dans le but de caractériser le paysage immunitaire du carcinome pulmonaire à petites cellules (LUAD). Une corrélation de Spearman a permis d'associer les gènes pronostiques candidats aux gènes de points de contrôle immunitaire. Afin d'évaluer les différences de réponse aux inhibiteurs de points de contrôle immunitaires (ICIs) chez les patients atteints de LUAD selon différents niveaux de risque, les données de score immunophénotypique (IPS) pour les traitements anti-PD-1 et anti-CTLA-4 ont été obtenues à partir de The Cancer Immunome Atlas (TCIA)18, et la base de données Tracking Tumor Immunophenotype (TIP)19 a été utilisée pour évaluer l'activité du cycle cancer-immunité, en comparant les scores correspondants entre les groupes à risque.
Analyse des mutations somatiques et de la sensibilité aux médicaments
L'outil de mutations TCGA a extrait les profils de mutations somatiques pour les cas TCGA-LUAD afin d'étudier la variation des motifs de mutation entre les groupes à risque. Le package maftools a permis le traitement et la visualisation des données de mutation. La charge mutationnelle tumorale (TMB) a été déterminée pour chaque échantillon et comparée entre les deux catégories de risque. Une analyse de sensibilité pharmacogénomique a été réalisée à l'aide du package pRRophetic selon la base de données Genomics of Drug Sensitivity in Cancer (GDSC)20. Les concentrations inhibitrices demi-maximales (IC50) pour les médicaments anticancéreux ont été prédites pour chaque patient atteint de LUAD, et les différences entre les groupes de risque ont été quantifiées à l'aide du test de Wilcoxon-Mann-Whitney.
Évaluation des niveaux d'expression de gènes pronostiques
Chaque jeu de données a permis d'évaluer les niveaux de transcrits de gènes candidats sélectionnés liés au résultat. Afin de relier l'expression génique au pronostic des patients, des seuils optimaux ont été déterminés à l'aide de la fonction surv_cutpoint du package R survminer. Sur la base de ces seuils, les cas de CPCA ont été classés en sous-groupes à expression élevée et faible pour les analyses de survie ultérieures.
De plus, cinq paires de tumeurs de carcinome pulmonaire à cellules adénoïdes (LUAD) et tissus normaux adjacents appariés ont été obtenues à l'hôpital thoracique d'Anhui, et une validation par RT-PCR quantitative a ensuite été réalisée. Chaque participant a fourni un consentement éclairé écrit. Six gènes pronostiques candidats (PDGFB, LDHA, ZEB2, FKBP4, DMD et S100B) ont été sélectionnés pour la validation par RT-PCR quantitative. L'ARN a été extrait à partir d'échantillons de tissus homogénéisés à l'aide d'un réactif d'extraction d'ARN, suivi d'une extraction au chloroforme et d'une précipitation à l'isopropanol. Un spectrophotomètre a permis la mesure de la concentration et de la pureté de l'ARN. La validation par RT-PCR quantitative des six gènes pronostiques candidats a été effectuée à l'aide d'un mélange maître PCR basé sur le SYBR Green sur un système de PCR en temps réel : une dénaturation initiale à 95 °C pendant 30 s, suivie de 40 cycles de 95 °C pendant 20 s, 55 °C pendant 20 s et 72 °C pendant 20 s. L'expression relative a été calculée et normalisée par rapport à la Glyceraldehyde-3-Phosphate Dehydrogenase (GAPDH) selon la méthode 2-ΔΔCt. Les détails de tous les réactifs et instruments sont fournis dans le Tableau des Matériaux.
Analyse statistique
Les analyses statistiques ont été réalisées à l'aide d'un logiciel de calcul et de représentation graphique statistiques. Le réseau d'interactions protéine-protéine a été visualisé à l'aide d'un logiciel d'analyse de réseaux. Après évaluation de la normalité, le test t de Student a été utilisé pour les variables continues normalement distribuées, et le test de Mann-Whitney U pour celles ne suivant pas une distribution normale.