Set di dati di benchmark FDA SEQC
Il set di dati di benchmark FDA SEQC (Gene Expression Omnibus accession number GSE47792) include dati di RNA-seq accoppiati generati utilizzando la piattaforma Illumina HiSeq 2000 con lunghezza di lettura di 100 nucleotidi7. Abbiamo utilizzato un sottoinsieme del set di dati SEQC-benchmark sequenziato in due siti: Beijing Genomics Institute (BGI) e Mayo Clinic (MAGGIO). Abbiamo usato quattro campioni (cioè A, B, C e D), ciascuno con quattro librerie di repliche preparate nei siti di sequenziamento. Il campione A contiene l’RNA di riferimento umano universale (UHRR), il campione B contiene l’RNA di riferimento del cervello umano (HBRR), il campione C contiene una miscela di A e B (75% A e 25% B) e il campione D contiene una miscela di A e B (25% A e 75% B). Abbiamo utilizzato i dati da due corsie di una singola cella di flusso per ogni replica del campione. Il SEQC ha anche fornito un set di dati di benchmark quantitativo PCR (qPCR) che include 20.801 geni analizzati con PrimePCR (Bio-Rad, Hercules, California). Ogni gene PrimePCR è stato analizzato una volta per ciascuno dei quattro campioni (cioè A, B, C e D). I set di dati e i campioni di benchmark FDA SEQC sono riassunti nelle tabelle supplementari S5 e S6.
Set di dati di neuroblastoma e adenocarcinoma polmonare
Abbiamo utilizzato un set di dati di neuroblastoma di 176 campioni (un sottoinsieme di un set di dati di 498 campioni più grande; indicato come SEQC-neuroblastoma in questo manoscritto) per valutare le prestazioni delle condotte RNA-seq in termini di previsione basata sull’espressione genica dell’esito della malattia. Questi campioni sono stati forniti dall’Ospedale pediatrico universitario di Colonia e sequenziati presso BGI utilizzando la piattaforma illumina48. Tutti i 176 campioni sono stati prelevati da pazienti ad alto rischio definiti come quelli con neuroblastoma di stadio 4 ed età > 18 mesi o con tumori amplificati da MYCN di qualsiasi stadio o età. Il set di dati SEQC-neuroblastoma è stato depositato al Gene Expression Omnibus con il numero di adesione GSE47792.
Abbiamo previsto due endpoint clinici: la sopravvivenza senza eventi (EFS), cioè il verificarsi di eventi come progresso, ricaduta o morte, e la sopravvivenza globale (OS), cioè la morte. Per entrambi gli endpoint, i pazienti sono stati suddivisi in due gruppi (ad alto rischio rispetto a basso rischio). I pazienti ad alto rischio hanno manifestato un evento o sono deceduti prima di una soglia di tempo di sopravvivenza predefinita, mentre i pazienti a basso rischio hanno manifestato un evento o sono deceduti dopo la soglia o il loro ultimo follow-up ha superato la soglia. Le soglie di tempo di sopravvivenza per EFS e OS erano rispettivamente di due e tre anni. Le soglie sono state scelte per bilanciare il numero di pazienti ad alto rischio e a basso rischio. I dettagli del set di dati SEQC-neuroblastoma sono forniti nella Tabella supplementare S9.
Abbiamo anche utilizzato un set di dati RNA-seq di adenocarcinoma polmonare di 87 campioni dal repository Cancer Genome Atlas (TCGA). L’endpoint di previsione era anche la sopravvivenza, e abbiamo utilizzato gli stessi criteri per definire gruppi ad alto rischio e a basso rischio con la soglia del tempo di sopravvivenza di due anni. La soglia di due anni è stata scelta per bilanciare il numero di pazienti ad alto rischio e basso rischio. I dettagli del set di dati TCGA-lung-adenocarcinoma sono forniti nella Tabella supplementare S10.
Filtraggio del set di dati del benchmark qPCR per produrre un set di riferimento di geni
A causa della variabilità nelle misurazioni del qPCR e dei disaccordi tra le piattaforme QPCR7, abbiamo filtrato il set di dati del benchmark qPCR per mantenere i geni che mostravano un comportamento “corretto”. Abbiamo quindi utilizzato questi geni per calcolare le metriche di benchmark (cioè accuratezza, precisione, affidabilità e riproducibilità). Tale processo di filtraggio è riassunto in Fig supplementare. S1.
Partendo dall’insieme iniziale di 20.801 geni analizzati con PrimePCR, abbiamo filtrato questi geni per mantenere solo geni che sono stati quantificati come diversi da zero (cioè rilevati) e con valori Ct (soglia di ciclo) ≤ 35 (35 indica il rilevamento di una sola molecola in un campione). Il filtraggio dei dati PrimePCR ha portato a 14.014 geni che hanno anche abbinato il trascrittoma AceView utilizzato per mappare il set di dati SEQC-benchmark RNA-seq.
Successivamente, abbiamo filtrato i geni 14,014 qPCR per mantenere solo i geni 12,610 che mostravano l’ordine di titolazione corretto (TO) e i rapporti di miscelazione attesi (EMR). I dettagli di questo processo sono nella sezione “Filtering qPCR genes by titation order and expected mixing ratios”.
Infine, poiché alcune metriche di benchmark come accuratezza e precisione sono sensibili a geni che esprimono zero o molto bassi, abbiamo ulteriormente selezionato geni che sono stati espressi come diversi da zero in tutte le repliche di tutti i campioni di tutti i siti di sequenziamento e tutte le pipeline 278 RNA – seq. Il set di riferimento finale contiene solo 10.222 geni qPCR (indicati come “tutti i geni”) che sono stati utilizzati per calcolare tutte e tre le metriche di riferimento per le pipeline RNA-seq.
Sulla base dello studio precedente, i geni con espressione inferiore hanno maggiori probabilità di essere incoerenti tra le pipeline49. Pertanto, abbiamo anche identificato un insieme di geni a bassa espressione nei geni 10,222 basati sull’espressione media qPCR dei campioni A, B, C e D. Il 20% più basso dei geni 10,222 (cioè, 2044 geni, indicati come “geni a bassa espressione”) sono stati utilizzati anche per calcolare lo stesso insieme di parametri di riferimento per le pipeline RNA-seq. Questo design ci ha permesso di studiare la capacità delle pipeline RNA-seq nella stima dell’espressione genica a bassa espressione.
Filtraggio dei geni qPCR per ordine di titolazione e rapporti di miscelazione attesi
I set di dati SEQC-benchmark (RNA-seq e qPCR) hanno proprietà uniche che consentono la valutazione della correttezza della quantificazione. Dopo aver identificato i geni qPCR rilevabili (cioè diversi da zero e Ct ≤ 35) e AceView-matched, abbiamo usato due metriche (TO ed EMR) per filtrare ulteriormente il set di dati qPCR di riferimento, lasciando solo i geni qPCR “corretti”. Le metriche TO e EMR acquisiscono proprietà di miscelazione uniche dei dati, ovvero
a Causa di questa proprietà, tutti i geni dovrebbero essere espressa in uno dei seguenti ordini, a seconda della relativa espressione di campioni A e B:
il set di qPCR geni che seguire la corretta titolazione ordine
Per un singolo replicare qPCR set di dati (ad es., il set di dati PrimePCR che abbiamo analizzato), la variabilità intrinseca di una singola misurazione qPCR può provocare alcuni geni falsi negativi che seguono il corretto ma non riescono ad essere identificati. Dal literature50, 51, il coefficiente di variazione per replicare le misurazioni qPCR è generalmente 15% o più grande, quindi abbiamo usato questo numero per regolare il margine per determinare se un gene segue la corretta A. Matematicamente, abbiamo calcolato l’intervallo di più e meno una deviazione standard da ogni misurazione qPCR e l’abbiamo usato come margine. Le equazioni rivedute per \({K}_{TO}\) sono le seguenti:
dove \(a=1.15, b=0.85\)
Oltre a, i campioni devono inoltre presentare un rapporto di miscelazione specifico. Dato che il rapporto tra i campioni A e B
EMR tra i campioni C e D è
e, infine, determina un insieme di geni che soddisfa le EMR criterio come segue:
dati di RNA-seq pipeline di analisi—mappatura, quantificazione, e la normalizzazione
Abbiamo studiato 278 RNA-seq condotte che comprendeva tredici sequenza di mappatura algorithms18,19,20,21,22,23,24,25,26,27,28,29, tre categorie di espressione quantificazione algorithms31,32,33, e sette espressione metodi di normalizzazione. Le tabelle supplementari S2-S4 riassumono tutte le opzioni considerate per ciascun componente della pipeline (mappatura delle sequenze, quantificazione delle espressioni e normalizzazione delle espressioni). Tredici mappatura algoritmi studiati sono Bowtie18, Bowtie219, BWA20, GSNAP21, Magic22 (un nuovo gasdotto sviluppato da NCBI per il SEQC progetto: ftp://ftp.ncbi.nlm.nih.gov/repository/acedb/Software/Magic), MapSplice23, Novoalign (un commercializzato pacchetto sviluppato da Novocraft: https://www.novocraft.com/products/novoalign/), OSA24, RUM25, STAR26, Subread27, TopHat28, e WHAM29. Alcuni usano la mappatura un-spliced delle letture al trascrittoma, e alcuni altri eseguono la mappatura spliced al genoma. Magic usa entrambi in parallelo e confronta la qualità di ogni allineamento per mantenere il meglio su più bersagli. Gli algoritmi di mappatura possono segnalare solo la mappatura univoca o consentire più posizioni di mappatura per lettura. Gli algoritmi di quantificazione includono metodi semplici basati sul conteggio (ad esempio, HTSeq31) e metodi probabilistici basati sulla distribuzione di Poisson applicati ai dati di mappatura genomica (ad esempio, Cufflinks32) o trascrittomica (ad esempio, RSEM33). La magia, il RUM e il sottotrame (cioè, featureCounts52) le pipeline includono metodi di quantificazione incorporati che rientrano nella categoria dei metodi semplici basati sul conteggio. I metodi di normalizzazione includono metodi di ridimensionamento semplici (cioè, frammenti per milione di frammenti mappati, frammenti per kilobase di lunghezza del gene per milione di frammenti mappati, mediana e quartile superiore), metodi di ridimensionamento robusti (cioè, espressione di log relativa e media tagliata di valori m) e metodi incorporati in pipeline specifiche (cioè, indice di espressione magica).
Sequence mapping
Abbiamo mappato le sequenze a ciascun riferimento in fasi successive utilizzando algoritmi di mappatura non-spliced o spliced. La mappatura un-spliced si riferisce ad algoritmi che allineano intere sequenze di lettura (ad esempio, Bowtie2, BWA e Novoalign) mentre la mappatura spliced si riferisce ad algoritmi che dividono le letture in segmenti per accogliere lunghi spazi vuoti o introni in una lettura (ad esempio, TopHat e MapSplice). Nella prima fase della mappatura un-spliced, abbiamo tentato di mappare tutte le sequenze accoppiate-end al riferimento ERCC / MT / rRNA (cioè, RNA esterno controlla sequenze Consorzio, il genoma mitocondriale, e sequenze di RNA ribosomiale). Tutte le coppie di lettura non mappate sono state quindi mappate al trascrittoma AceView. Infine, tutte le coppie di lettura che non sono state mappate ai riferimenti ERCC/MT/rRNA o AceView sono state mappate al riferimento del genoma umano. Le coordinate di mappatura trascrittomica sono state poi tradotte in coordinate di mappatura genomica e fuse con i risultati di mappatura del genoma umano per produrre i risultati finali (Fig. S21, pannello sinistro). Abbiamo usato Bowtie2 come mappatore per il primo passo di tutte le pipeline di mappatura giuntate (Fig. S21, pannello di destra). Algoritmi di mappatura impiombato direttamente mappato legge al genoma umano (ad esempio, MapSplice e GSNAP) o mappato intero non-impiombato legge al trascrittoma e poi unito questi risultati di mappatura con i risultati di mappatura impiombato del restante legge al genoma umano (ad esempio, TopHat e T). La tabella supplementare S2 riassume tutti gli strumenti di mappatura studiati in questo studio.
Bowtie2, GSNAP, Novoalign, TopHat e WHAM consentono il controllo del numero di mappature riportate per coppia di lettura. Per impostazione predefinita, questi algoritmi in genere riportano una singola posizione di mappatura migliore per coppia di lettura. Tuttavia, alcuni algoritmi di quantificazione possono utilizzare informazioni su più posizioni di mappatura ambigue per migliorare la stima dell’espressione genica. Così, oltre al reporting single-hit, abbiamo anche generato risultati di mappatura che hanno riportato fino a 200 colpi per lettura (multi-hit). Abbiamo anche incluso la pipeline di mappatura Bowtie con parametri di mappatura specifici per la quantificazione con RSEM, come descritto nella seguente sezione33.
Le opzioni della riga di comando per tutti gli strumenti di allineamento delle sequenze sono dettagliate nella Nota supplementare 1.
Quantificazione dell’espressione genica
La fase di quantificazione comprendeva tre categorie di quantificatori-quantificatori basati sul conteggio (cioè, HTSeq e quantificatori incorporati per le pipeline Magic, RUM e Subread), quantificatori basati su modelli di probabilità per la mappatura genomica (cioè, Gemelli), e quantificatori basati su modelli di probabilità per la mappatura trascrittomica (cioè, RSEM). Le caratteristiche principali di questi quantificatori sono riassunte nella tabella supplementare S3. Gemelli è un quantificatore basato su modello di Poisson che stima le probabilità di assegnazione di lettura in base alle informazioni di allineamento32. È in grado di assemblare trascritti e quantificare espressioni geniche o trascrittive. In questo studio, abbiamo disabilitato la funzione di assemblaggio e fornito il file genome annotation GTF come riferimento di quantificazione. HTSeq è un quantificatore naïve basato sul conteggio che assegna letture mappate a genes31. HTSeq è in grado di quantificare l’espressione genica, ma non di trascrivere l’espressione. RSEM è anche un quantificatore basato sul modello di Poisson che è simile nel concetto a Cufflinks33. Le informazioni da letture multi-hit sono importanti sia per i gemelli che per RSEM. Questi algoritmi utilizzano multi-hit leggere le informazioni per stimare più accuratamente gene o espressione trascrizione.
I risultati della mappatura delle pipeline di allineamento non erano sempre compatibili con le tre categorie di quantificatori. Gemelli richiede che i risultati di allineamento sono ordinati per coordinate di allineamento e letture multi-hit sono contrassegnati con il tag ‘ NH ‘ nel campo attributo del file SAM. HTSeq richiede che i risultati dell’allineamento siano ordinati per nomi letti e che il tag’ NH ‘ sia assente nel file SAM. RSEM quantifica solo la mappatura trascrittomica, cioè legge mappata e riportata in coordinate trascrittomiche. Inoltre, RSEM gestisce solo allineamenti non gappati. Pertanto, è necessario filtrare per rimuovere gli allineamenti gappati. A causa di questi requisiti, abbiamo pre-elaborato tutti i risultati di allineamento prima della quantificazione. In sintesi, venti pipeline di allineamento, tra cui pipeline spliced, un-spliced, single-hit e multi-hit, erano adatte per la quantificazione basata sul conteggio. Sedici condotte di allineamento erano adatte per gemelli e solo dieci erano adatte per RSEM. RSEM è specificamente progettato per funzionare bene con Bowtie. Pertanto, abbiamo incluso anche questa pipeline di mappatura e quantificazione incorporata.
Le opzioni della riga di comando per tutti gli strumenti di quantificazione sono dettagliate nella Nota supplementare 1.
Normalizzazione dell’espressione genica
La normalizzazione dei dati RNA-seq consente il confronto tra campioni. Generalmente, i metodi di normalizzazione correggono la dimensione della libreria (cioè il numero totale di letture in un campione), che è la fonte primaria della varianza tra campioni. Abbiamo studiato sette metodi di normalizzazione-frammenti per milione di frammenti mappati (FPM), frammenti per kilobase di lunghezza del gene per milione di frammenti mappati (FPKM), mediana (Med.), quartile superiore (UQ), espressione di log relativa (RLE), media tagliata di M-values (TMM) e indice di espressione (EIndex, che è specifico per la pipeline magica) (vedere Tabella supplementare S4). Descriviamo ciascuno di questi metodi di normalizzazione in base alla seguente descrizione matematica del set di dati SEQC-benchmark.
abbiamo definito l’insieme dei geni presenti per essere
e il finale attuale gene set
Abbiamo usato lo stesso insieme di gens presenti per tutti i metodi di normalizzazione per una pipeline RNA-seq.
Il numero totale di geni presenti in un dato campione s e replicare n è
e la media del numero totale di geni presenti per tutti i dati da un singolo sito è
Così, abbiamo definito FPM-espressione normalizzata per ogni campione s, replicare n, e gene k come
Mediana e superiore quartile normalizzato di espressione per ogni campione s, replicare n, e gene k sono quindi definiti come
Per la normalizzazione FPKM, abbiamo definito la lunghezza di un gene k come \(\ell_{k}\), che è la lunghezza dell’unione di tutti gli esoni relativi al gene come definito dal trascrittoma AceView. La formulazione originale di FPKM arbitrariamente utilizzato fattori di scala di 1 × 103 per la lunghezza del gene e 1 × 106 per il numero totale di frammenti mappati. Al fine di mantenere la gamma dinamica comparabile tra tutti i metodi di normalizzazione, abbiamo invece scalato dalla lunghezza media del gene e dal conteggio totale medio per tutti i geni presenti. La lunghezza media di tutti i geni presenti è
Così, riproporzionato FPKM normalizzato di espressione per ogni campione s, replicare n, e gene k è
I metodi di normalizzazione TMM e RLE sono simili alla normalizzazione FPM, ma introducono un fattore di ridimensionamento aggiuntivo per regolare le dimensioni della libreria. Abbiamo usato il pacchetto edgeR in R per stimare un fattore di scala per ogni replica del campione36, 53. Il metodo TMM seleziona una libreria di riferimento da un pool di librerie di replica campione e quindi calcola i rapporti di espressione di log gene-wise (valori M) e i valori di espressione di log media gene-wise (valori A) tra la libreria di destinazione e la libreria di riferimento. I numeri estremi nei valori M e A vengono tagliati e il fattore di scala per la libreria di destinazione è la media ponderata dei valori M rimanenti. Il metodo RLE determina un fattore di scala definendo innanzitutto la libreria mediana come media geometrica genica attraverso i replicati campioni35. Il rapporto mediano tra ciascuna libreria di destinazione e la libreria mediana viene preso come fattore di scala. TMM – e RLE normalizzato di espressione per ogni campione s, replicare n, e gene k è definito come:
dove \(\hat{f}_{s,n}^{TMM}\) e \(\hat{f}_{s,n}^{RLE}\) sono il fattore di scala, per esempio, di s, di replicare n.
RNA-seq pipeline di metriche di performance
Benchmark di metriche per la RNA-seq oleodotti, sono riassunti nella Tabella Supplementare S7.
Precisione misurata come deviazione da qPCR riferimenti
Dato che i campioni A e B, l’assoluta log-rapporto di deviazione di RNA-seq-based espressione da qPCR a base di espressione di un gene k è
e la metrica di precisione finale è stata definita come la mediana di tutti \({\Delta }_{{\frac{A}{B},k}}\), \(k = 1 \ldots K\).
Precisione misurata come variazione dell’espressione genica attraverso le librerie replicate
Abbiamo calcolato il coefficiente di variazione (CoV) per ciascun gene e ciascun campione attraverso quattro librerie replicate come segue:
Affidabilità misurata come intra-esempio di correlazione di espressione genica
L’affidabilità di un sistema di misura può essere valutata con il coefficiente di correlazione intraclasse (ICC)54,55. IC è applicabile alle misurazioni che possono essere organizzate in gruppi e descrive come le misurazioni simili dello stesso gruppo sono tra loro. La moderna definizione IC prende in prestito il quadro di analisi della varianza (ANOVA), o più specificamente ANOVA con effetti casuali55. Il tipo di ANOVA dipende dal progetto sperimentale e generalmente segue la definizione nell’articolo di Shrout pubblicato nel 197955. IC(1,1) e IC (1,k) si basano sul modello di effetti casuali unidirezionali e sono applicabili al caso in cui ciascun gruppo sia valutato da un diverso insieme di valutatori k selezionati casualmente da una popolazione più ampia di valutatori. IC(2,1) e IC (2,k) si basano sul modello di effetti casuali a due vie e sono applicabili al caso in cui un campione casuale di k raters è preselezionato da una popolazione più ampia e ogni rater valuta ogni gruppo esattamente una volta (cioè, ogni rater valuta n gruppi complessivamente). IC(3,1) eC (3,k) si basano sul modello a effetti misti a due vie e sono applicabili al caso in cui ciascun gruppo sia valutato da ciascuno degli stessi valutatori k, che sono gli unici valutatori della popolazione. Il secondo parametro in IC (,) indica se l’IC è quello di misurare l’affidabilità di una singola misura o la media di k misurazioni.
Per il set di dati di benchmark SEQC con librerie replicate per ciascun campione, IC(1,1) o IC(1,k) hanno fissato il nostro obiettivo poiché, per un gene g, l’espressione genica di librerie replicate per campioni diversi (o gruppi diversi nel contesto precedente) non è stata valutata esattamente nelle stesse condizioni (o valutata dagli stessi valutatori nel contesto precedente). Abbiamo scelto di utilizzare IC (1, k) come librerie di replica sono disponibili per la maggior parte degli esperimenti. Matematicamente, un one-way modello a effetti casuali possono essere formulati come
Abbiamo calcolato ICC per ogni gene k, \(k = 1 \ldots K\), e quindi utilizzata la mediana di tutti ICCs come misura finale di affidabilità.
Abbiamo anche studiato altre metriche potenziali, come la riproducibilità, che è definita come la correlazione di Spearman tra due librerie replicate dello stesso campione (Nota supplementare 2). La correlazione Spearman variava da 0.da 993 a 0,996 (Fig. S8) utilizzando AllGenes. Abbiamo scartato la metrica di riproducibilità a causa della gamma dinamica relativamente piccola.
Valutando l’utilità delle metriche di benchmark per la selezione della pipeline RNA-Seq
Abbiamo classificato le pipeline RNA-seq in base al rango medio delle tre metriche di benchmark (cioè accuratezza, precisione e affidabilità). Abbiamo quindi valutato l’utilità del benchmark metriche esaminando se di buona esecuzione e poveri-l’esecuzione di tubazioni identificato sulla base di benchmark di metriche stati informativo per la deduzione delle prestazioni di espressione genica di predizione dell’esito della malattia e la significatività statistica del paziente stratificazione per tutti gli endpoint clinici (cioè, il SEQC-neuroblastoma EFS e OS endpoint e il TCGA-polmone-adenocarcinoma endpoint di sopravvivenza).
In primo luogo, per le 278 pipeline rappresentative RNA-seq applicate al set di dati SEQC-benchmark, abbiamo calcolato il rank medio utilizzando un sottoinsieme delle metriche di benchmark come indicatore di performance finale per ogni pipeline. In totale, abbiamo avuto 6 metriche (3 benchmark metriche × 2 set di gene ), e abbiamo studiato 12 sottoinsiemi (4 × 3) del 6 metriche utilizzando i seguenti criteri:
- (1)
Quattro combinazioni di tre benchmark metriche con almeno due in un sottoinsieme—una combinazione di tutti e tre i benchmark metriche, tre combinazioni con due delle tre benchmark metriche.
- (2)
Tre sottoinsiemi formati da metriche derivate da tutti i geni, quelli derivati da geni a bassa espressione o da una combinazione di entrambi.
In secondo luogo, per ciascuna delle 278 pipeline rappresentative RNA-seq (156 per l’endpoint di sopravvivenza TCGA-lung-adenocarcinoma), abbiamo calcolato AUC e MCC di convalida incrociata nidificate, come descritto nella sezione “Metodo” “Neuroblastoma and lung adenocarcinoma predictive modeling”, risultando in 834 (468 per l’endpoint di sopravvivenza TCGA-lung-adenocarcinoma) valori di AUC e MCC per ciascun endpoint clinico (cioè, 278 pipeline × 3 classificatori, o 156 pipeline × 3 classificatori) (Tabelle supplementari S11, S12). Abbiamo anche modellato le funzioni di sopravvivenza usando l’analisi Kaplan–Meier per ogni pipeline, come descritto nella sezione “Metodo” “Analisi di sopravvivenza Kaplan–Meier”. Per ogni pipeline RNA-seq, abbiamo riassunto le prestazioni della previsione basata sull’espressione genica dell’esito della malattia utilizzando sia l’AUC media che l’MCC tra i classificatori sia il tasso di successo della stratificazione del paziente (cioè, separazione statisticamente significativa di due curve Kaplan-Meier) in tutte le iterazioni e classificatori nel framework di convalida incrociata nidificato.
Infine, abbiamo identificato le principali pipeline con buone prestazioni del 10% e le pipeline con cattive prestazioni del 10% in base al rango medio di un sottoinsieme delle tre metriche di benchmark. La prestazione di previsione corrispondente (cioè, AUC e MCC) delle pipeline di buone prestazioni è stato testato contro quello delle pipeline di scarso rendimento utilizzando il test di Wilcoxon rank-sum unilaterale con l’ipotesi nulla che la mediana del primo gruppo non fosse più grande di quella del secondo gruppo.
Modellazione predittiva del neuroblastoma e dell’adenocarcinoma polmonare
Abbiamo valutato le prestazioni di 278 pipeline RNA-seq in termini di processo decisionale basato sull’espressione genica utilizzando i dati SEQC-neuroblastoma 48. Il set di dati SEQC-neuroblastoma e gli endpoint clinici associati sono riassunti nella Tabella supplementare S9. Le pipeline RNA-seq sono state valutate in termini di previsione degli esiti dei pazienti con neuroblastoma per due endpoint clinici utilizzando validazione incrociata nidificata (Fig. S13) 56,57. Allo stesso modo abbiamo anche valutato le prestazioni di 156 pipeline RNA-seq applicate al set di dati TCGA-lung-adenocarcinoma per predire l’esito della malattia. Il set di dati TCGA-lung-adenocarcinoma e l’endpoint clinico associato sono riassunti nella Tabella supplementare S10.
La convalida incrociata nidificata comporta la formazione e il test di un modello di previsione ottimale. Ciò si ottiene utilizzando la triplice ottimizzazione o validazione incrociata interna, applicata al sottoinsieme di allenamento dalla validazione incrociata esterna di cinque volte. Una volta identificati i parametri finali del modello di previsione ottimale (ad esempio, gli iperparametri del classificatore e la dimensione della caratteristica), il modello finale viene addestrato utilizzando l’intero sottoinsieme di allenamento e quindi testato utilizzando la piega rimanente dalla convalida incrociata esterna di cinque volte. Questo processo è stato ripetuto per dieci iterazioni. Abbiamo condotto la convalida incrociata nidificata separatamente per ciascuno dei tre classificatori (cioè, adaptive boosting, logistic regression e support vector machines) e ha utilizzato il metodo di selezione delle funzionalità Minimum redundancy, Maximum Relevance (mRMR) per scegliere le dimensioni ottimali delle funzionalità tra 5 e 40 con la dimensione del passo di 558.
Kaplan–Meier survival analysis
Per ogni pipeline e classificatore RNA-seq (cioè, 278 pipeline × 3 classificatori per gli endpoint SEQC-neuroblastoma e 156 pipeline × 3 classificatori per l’endpoint di sopravvivenza TCGA-lung-adenocarcinoma), abbiamo modellato le funzioni di sopravvivenza di Kaplan–Meier in base alle etichette previste di ciascun campione. Abbiamo quindi utilizzato il test log-rank a due code per determinare se le curve di sopravvivenza stimate per ciascun gruppo di pazienti previsto fossero statisticamente diverse.
l’Analisi della varianza e il calcolo del contributo di ogni RNA-seq pipeline fattore generale pipeline varianza
Abbiamo utilizzato l’analisi della varianza (ANOVA) per determinare se ogni RNA-seq pipeline fattore che contribuisce in modo significativo alla varianza di ciascuna delle tre benchmark metriche (cioè, l’accuratezza, la precisione e l’affidabilità), come pure per la varianza di stima delle prestazioni (cioè, l’AUC e la MCC). Per ciascuna delle tre metriche di benchmark, abbiamo utilizzato un modello lineare (funzione R “lm”) per adattare i dati di tutte le 278 pipeline utilizzando la metrica come variabile dipendente e i fattori della pipeline RNA-seq come variabili categoriali indipendenti. Abbiamo considerato i seguenti fattori come variabili categoriche indipendenti-algoritmo di mappatura, strategia di mappatura (cioè, spliced vs. un-spliced), reporting di mappatura (cioè, single-hit vs. multi-hit), algoritmo di quantificazione e algoritmo di normalizzazione. Abbiamo incluso tutti i fattori e le loro interazioni bidirezionali nel modello lineare. Per ciascuno degli endpoint di previsione, abbiamo applicato la stessa tecnica per adattare i dati di tutte le 278 pipeline utilizzando AUC media o MCC come variabile dipendente e lo stesso insieme di fattori di pipeline RNA-seq come variabili categoriali indipendenti. Abbiamo quindi condotto l’ANOVA sul modello lineare (funzione R “anova”). ANOVA calcola una “somma di quadrati” (cioè varianza) attribuita a ciascun fattore o interazione e utilizza un test F per determinare se la varianza è statisticamente significativa. Abbiamo calcolato la percentuale che ogni fattore o interazione contribuisce alla varianza totale calcolando il rapporto tra “somma dei quadrati” per ciascun fattore e la somma totale dei quadrati.
Analisi di regressione
Abbiamo studiato la relazione tra profili di allineamento o caratteristiche di distribuzione dell’espressione genica e metriche di benchmark. I profili di allineamento includevano il numero totale di frammenti mappati, il numero totale di letture che coprono la regione intronica, il numero totale di letture con inserimenti o eliminazioni, il numero totale di letture perfettamente corrispondenti, il numero totale di letture con al massimo una mancata corrispondenza e il numero di disallineamenti per lettura mappata. Ogni algoritmo di allineamento era rappresentato dalle statistiche medie su 2 siti di sequenziamento, 4 campioni, 4 librerie di repliche e 2 corsie. Utilizzando il pacchetto “MASS” in R, abbiamo adottato l’approccio M-estimation con ponderazione Huber per adattare modelli di regressione lineare robusti tra una variabile dipendente (prestazioni metriche di riferimento) e una variabile esplicativa (un profilo di allineamento). La stima M con l’approccio di ponderazione Huber è un metodo di regressione che è robusto in presenza di valori anomali. Le caratteristiche di distribuzione dell’espressione genica includevano il quartile inferiore, il quartile medio, il quartile superiore, il massimo, l’intervallo interquartile, la deviazione standard, l’asimmetria, la curtosi e l’entropia di una distribuzione di espressione genica. Abbiamo utilizzato la stessa stima M con l’approccio di ponderazione Huber per adattarsi a modelli di regressione lineare robusti e quindi abbiamo riportato l’errore standard residuo per ciascun modello.
Disclaimer
Le opinioni presentate in questo articolo non riflettono necessariamente l’opinione o la politica attuale o futura della Food and Drug Administration statunitense. Qualsiasi menzione di prodotti commerciali è per chiarimenti e non inteso come approvazione.