$$\rightleftharpoonup{xx}$$
$$\longleftharp{xx}$$,
$$\longrightharp{xx}$$,
Questo studio ha utilizzato solo dataset pubblicamente disponibili e de-identificati e non ha coinvolto sperimentazioni dirette su uomini o animali; Pertanto, non era necessaria ulteriore approvazione del comitato etico e un consenso informato.
Download e elaborazione dei dati
I dati di sequenziamento dell'RNA e le relative informazioni cliniche per il carcinoma a cellule squamose polmonari (LUSC) sono stati ottenuti dal database The Cancer Genome Atlas (TCGA) tramite il portale dati Genomic Data Commons nell'ambito del progetto TCGA-LUSC. La matrice di espressione TCGA-LUSC utilizzata in questo studio si basava sui valori FPKM. I valori di espressione genica sono stati trasformati e normalizzati prima delle analisi a valle. Le variabili cliniche includevano età, sesso, stadio tumorale, stadio patologico della TNM, grado, tempo di sopravvivenza e stato di sopravvivenza quando disponibili. Un totale di 489 casi TCGA-LUSC sono stati inizialmente prelevati, e 381 pazienti con profili di espressione completi e informazioni di sopravvivenza complessiva sono stati inclusi nella costruzione del modello prognostico e nella valutazione interna.
Dataset indipendenti di validazione sono stati scaricati dal database Gene Expression Omnibus (GEO) il 3 gennaio 2026. GSE30219 si basava sulla piattaforma GPL570 e includeva 307 pazienti LUSC con informazioni disponibili sull'espressione e sulla sopravvivenza complessiva. GSE37745 si basava anche sulla piattaforma GPL570 e includeva 196 pazienti LUSC con informazioni disponibili sull'espressione e sulla sopravvivenza complessiva. GSE57148 si basava sulla piattaforma GPL11154 e includeva 91 tessuti polmonari normali e 98 tessuti polmonari provenienti da pazienti con broncopatia cronica ostruttiva (BPCO), per un totale di 189 campioni. GSE57148 veniva utilizzato per identificare geni differenzialmente espressi associati alla BPCO, mentre GSE30219 e GSE37745 venivano usati come coorti di validazione esterna indipendenti.
Per i dataset GEO, l'annotazione delle sonde veniva eseguita utilizzando il corrispondente pacchetto di annotazione R e i file di annotazione della piattaforma. Gli identificatori delle sonde venivano convertiti in simboli genici ufficiali. Quando più sonde erano mappate sullo stesso gene, la sonda con il valore medio di espressione più alto veniva mantenuta per rappresentare quel gene. Non sono stati rilevati geni del modello mancanti nei dataset di validazione dopo l'abbinamento gene-simbolo. TCGA-LUSC, GSE30219 e GSE37745 sono stati analizzati come coorti generali di LUSC perché lo stato di comorbidità a livello di paziente non è stato confermato nelle annotazioni utilizzate per questa analisi.
Poiché i dataset TCGA e GEO venivano generati utilizzando diverse piattaforme di espressione, la normalizzazione cross-platform e la correzione degli effetti batch venivano eseguite utilizzando metodi standard di preprocessing basati su R prima dell'applicazione del modello. Il modello di rischio è stato addestrato nella coorte TCGA-LUSC e poi valutato indipendentemente in ciascuna coorte esterna GEO, invece di fondere direttamente tutte le coorti. L'analisi differenziale dell'espressione veniva eseguita utilizzando il pacchetto limma R. I geni differenzialmente espressi a BPCO in GSE57148 sono stati screenati utilizzando |log2FC| > 0,263 e P < 0,05. Questa soglia log2FC corrisponde a un cambiamento di circa 1,2 volte ed è stata utilizzata come criterio di screening esplorativo per mantenere geni potenzialmente rilevanti associati alla PANoptosi. Un totale di 277 geni associati alla PANoptosi sono stati selezionati da studi precedentemente pubblicati su apoptosi - - - pirottosi, necroptosi - e PANoptosi e sono forniti nella Tabella Supplementare 1.
Analisi dell'arricchimento funzionale dei geni
Per chiarire le implicazioni funzionali dei geni PANoptosi associati alla COPD selezionati, sono state effettuate analisi di arricchimento Gene Ontology (GO) e Kyoto Encyclopedia of Genes and Genomes (KEGG) utilizzando il pacchetto clusterProfiler R e l'org. Hs.eg.db pacchetto di annotazione. Sono state valutate le categorie di processo biologico GO, componente cellulare e funzione molecolare. L'analisi di arricchimento delle vie KEGG è stata condotta per identificare le vie di segnalazione associate ai geni selezionati. I valori P sono stati regolati per test multipli utilizzando il metodo del tasso di falsa scoperta di Benjamini-Hochberg, dove applicabile. I termini di arricchimento con P < 0,05 sono stati considerati statisticamente significativi in questa analisi esplorativa. I plot di arricchimento sono stati generati usando ggplot2.
Analisi di clustering non supervisionata dei modelli di espressione genica associati alla PANoptosi
Per esplorare l'eterogeneità molecolare associata all'espressione genica associata alla PANoptosi in LUSC, è stato effettuato un clustering di consenso utilizzando il pacchetto Consensus ClusterPlus R. I campioni di LUSC sono stati raggruppati secondo i profili di espressione dei geni PANoptosi associati alla sopravvivenza. Il clustering gerarchico è stato applicato con la distanza di correlazione di Pearson. Il numero massimo di cluster è stato impostato a sei e sono state effettuate 1.000 iterazioni di ricampionamento per valutare la robustezza del clustering. Il numero ottimale di cluster è stato determinato valutando la matrice di consenso, la curva della funzione di distribuzione cumulativa, il grafico dell'area delta e l'interpretabilità biologica dei gruppi risultanti. Sulla base di questi criteri, k = 2 è stato selezionato per l'analisi a valle. Le differenze di sopravvivenza tra i due gruppi molecolari sono state valutate utilizzando l'analisi di Kaplan-Meier e il test log-rank.
Analisi delle differenze del microambiente immunitario tra sottotipi
Per confrontare le caratteristiche del microambiente immunitario tra sottotipi molecolari, l'infiltrazione delle cellule immunitarie è stata stimata utilizzando l'algoritmo di deconvoluzione CIBERSORT con la matrice di firma leucocitaria LM22. L'analisi è stata eseguita in R utilizzando i pacchetti e1071 e preprocessCore. I valori P della permutazione CIBERSORT sono stati registrati per valutare l'affidabilità delle stime di deconvoluzione. Poiché questo studio era esplorativo e basato su dati trascricromici retrospettivi, le differenze tra cellule immunitarie sono state interpretate come modelli di infiltrazione immunitaria dedotti computazionalmente piuttosto che come misurazioni cellulari dirette.
L'algoritmo ESTIMATE è stato utilizzato per calcolare il punteggio stromale, il punteggio immunitario, il punteggio ESTIMATE e la purezza tumorale per ciascun campione tumorale. Il GSVA è stato applicato per stimare i punteggi di arricchimento a livello di percorso basati su set genetici selezionati. Le differenze a gruppo nelle frazioni delle cellule immunitarie, nei geni dei punti di controllo immunitario, nei geni della famiglia HLA e nei punteggi derivati da ESTIMATE sono state valutate utilizzando test non parametrici. Per confronti immunitari multipli, è stata applicata la correzione di Benjamini-Hochberg dove appropriato; le analisi riportate utilizzando valori nominali di P sono state interpretate come esplorative. L'analisi di correlazione di rango di Spearman è stata utilizzata per valutare le associazioni tra espressione genica e marcatori immunocorrelati, con sia i coefficienti di correlazione che i valori P riportati quando applicabile. Le differenze tra i trascritori HLA sono state interpretate come cambiamenti trascrizionali correlati alla presentazione dell'antigene piuttosto che come evidenza funzionale diretta di una maggiore capacità di presentazione dell'antigene.
Stabilimento di una firma prognostica correlata ai geni associati alla PANoptosi
La coorte TCGA-LUSC con profili di espressione completi e informazioni sulla sopravvivenza complessiva è stata utilizzata per la costruzione prognostica del modello. Tra i 489 casi TCGA-LUSC inizialmente recuperati, 381 pazienti con dati completi sulla sopravvivenza complessiva sono stati inclusi nell'analisi prognostica. Questi pazienti sono stati divisi casualmente in una coorte di addestramento e una coorte interna di test con un rapporto di 7:3. La randomizzazione stratificata è stata eseguita in base allo stato di sopravvivenza per mantenere una distribuzione comparabile degli eventi di sopravvivenza tra le coorti di addestramento e test.
Nella coorte di addestramento, la regressione dei rischi proporzionali di Cox univariata è stata utilizzata per la prima volta per valutare l'associazione tra ciascun gene associato alla PANoptosi candidato e la sopravvivenza complessiva. I geni con P < 0,05 erano considerati geni prognostici candidati e successivamente inseriti nella regressione LASSO Cox utilizzando il pacchetto glmnet R. È stata utilizzata una validazione incrociata dieci volte per selezionare il parametro di penalizzazione ottimale e ridurre l'overfitting. Sulla base dei coefficienti di regressione LASSO Cox e dei corrispondenti valori di espressione genica normalizzati, è stato calcolato un punteggio di rischio individualizzato per ciascun paziente utilizzando la formula:
Punteggio di rischio = Σ(coefi × Xi)
dove coefi rappresenta il coefficiente di regressione di ciascun gene selezionato e Xi rappresenta il valore di espressione normalizzato del gene corrispondente. Il modello prognostico finale conteneva 12 geni: GSDMD, IL18, NFKBIA, PIK3CA, IL1B, BIRC3, MCL1, PSMB10, LMNA, CFLAR, IL1R1 e AKT3. L'equazione completa basata sul punteggio di rischio basata sui coefficienti è fornita nella Tabella Supplementare 2.
Il punteggio mediano di rischio nella coorte di formazione è stato utilizzato come cutoff per classificare i pazienti in gruppi ad alto rischio e a basso rischio. La stessa formula di punteggio di rischio è stata applicata alla coorte interna di test e alle coorti di validazione esterne GSE30219 e GSE37745. Sono state utilizzate analisi di sopravvivenza Kaplan-Meier, test log-rank e analisi della curva funzionante del ricevitore dipendenti dal tempo per valutare le prestazioni del modello. Poiché i dataset di validazione sono stati generati utilizzando piattaforme microarray e non contenevano annotazioni di comorbidità COPD confermate, la validazione esterna è stata interpretata come valutazione retrospettiva in coorti LUSC indipendenti piuttosto che come validazione in pazienti clinicamente confermati con LUSC COPD-COD.
Analisi della previsione della sensibilità ai farmaci
La sensibilità al farmaco è stata stimata utilizzando il package pRRophetic R, che prevede la risposta al farmaco a partire dai profili di espressione genica tumorale basandosi sui dati di riferimento farmacogenomici del database Genomics of Drug Sensitivity in Cancer. Per ogni campione di paziente sono stati calcolati valori previsti di concentrazione inibitoria a metà massima. Le matrici di espressione venivano elaborate secondo i requisiti di input pRRophetic, e la correzione batch-effect veniva eseguita utilizzando il flusso di lavoro standard compatibile con pRRophetic. I valori previsti di IC50 sono stati riportati sulla scala di output pRRophetica.
Otto agenti candidati, tra cui sorafenib, gefitinib, bleomicina, bosutinib, etoposide, lenalidomide, camptothecina e metotrexato, sono stati valutati come pannello esplorativo di sensibilità ai farmaci. Le differenze nei valori previsti di IC50 tra gruppi ad alto e basso rischio sono state confrontate utilizzando il test di somma dei rangi di Wilcoxon. Questi risultati sono stati interpretati come stime computazionali della sensibilità ai farmaci piuttosto che come risposta chemoterapica clinica misurata o resistenza ai farmaci confermata sperimentalmente.
Analisi statistica
Tutte le analisi statistiche sono state condotte utilizzando il software R. Le variabili continue tra due gruppi sono state confrontate utilizzando il test di Wilcoxon a somma dei rangi, mentre confronti tra più di due gruppi sono stati effettuati utilizzando il test di Kruskal-Wallis quando appropriato. La sopravvivenza complessiva è stata definita come il principale punto di riferimento della sopravvivenza. Sono state generate curve di sopravvivenza di Kaplan-Meier per confrontare le differenze di sopravvivenza tra i gruppi, e la significatività statistica è stata valutata utilizzando il test log-rank. Univariata e multivariata
Sono state utilizzate analisi di regressione dei rischi proporzionali di Cox per valutare le associazioni prognostiche tra variabili cliniche, gruppo di rischio e sopravvivenza complessiva. Sono state considerate variabili con rilevanza clinica o significatività statistica nell'analisi di Cox univariata per la regressione di Cox multivariata. L'assunzione dei rischi proporzionali è stata valutata utilizzando i residui di Schoenfeld. Le variabili cliniche mancanti sono state gestite utilizzando analisi completa del caso per la regressione di Cox e la costruzione del nomogramma. La collinearità tra le variabili cliniche è stata valutata prima della modellazione multivariata.
Curve ROC dipendenti dal tempo sono state utilizzate per valutare la performance predittiva del modello di rischio per la sopravvivenza complessiva su 1, 3 e 5 anni. Un nomogramma è stato costruito utilizzando variabili mantenute nel modello multivariato o variabili con sufficiente disponibilità clinica. Sono stati utilizzati grafici di calibrazione per confrontare le probabilità di sopravvivenza complessiva previste e osservate. L'analisi della curva decisionale è stata effettuata come valutazione esplorativa del potenziale beneficio netto attraverso probabilità di soglia selezionate.
L'analisi di correlazione dei ranghi di Spearman è stata utilizzata per valutare le associazioni tra espressione genica e caratteristiche immunitarie. Coefficienti di correlazione e valori P sono stati riportati dove applicabile. Per confronti multipli, è stata applicata la correzione del tasso di falsa scoperta di Benjamini-Hochberg dove appropriato. Le analisi riportate utilizzando valori nominali di P sono state considerate esplorative. Un valore P bilaterale < 0,05 è stato considerato statisticamente significativo.