Impacto de los algoritmos de análisis de datos ARN-seq en la estimación de la expresión génica y la predicción descendiente

Conjuntos de datos de referencia FDA SEQC

El conjunto de datos de referencia FDA SEQC (Número de acceso general de expresión génica GSE47792) incluye datos de ARN-seq de extremo emparejado generados utilizando la plataforma Illumina HiSeq 2000 con una longitud de lectura de 100 nucleótidos 7. Utilizamos un subconjunto del conjunto de datos de referencia SEQC secuenciado en dos sitios: el Instituto de Genómica de Beijing (BGI) y la Clínica Mayo (MAY). Se utilizaron cuatro muestras (es decir, A, B, C y D), cada una con cuatro bibliotecas replicadas preparadas en los sitios de secuenciación. La muestra A contiene el ARN Universal de Referencia Humana (RDH), la muestra B contiene el ARN de Referencia del Cerebro Humano (RBRH), la muestra C contiene una mezcla de A y B (75% A y 25% B), y la muestra D contiene una mezcla de A y B (25% A y 75% B). Utilizamos datos de dos carriles de una sola celda de flujo para cada réplica de muestra. El SEQC también proporcionó un conjunto de datos de referencia de PCR cuantitativa (qPCR) que incluye 20.801 genes analizados con PrimePCR (Bio-Rad, Hercules, California). Cada gen PrimePCR se analizó una vez para cada una de las cuatro muestras (es decir, A, B, C y D). Los conjuntos de datos y muestras de referencia de la FDA SEQC se resumen en las Tablas complementarias S5 y S6.

Conjuntos de datos de neuroblastoma y adenocarcinoma de pulmón

Utilizamos un conjunto de datos de neuroblastoma de 176 muestras (un subconjunto de un conjunto de datos de 498 muestras más grande; se denomina SEQC-neuroblastoma en este manuscrito) para evaluar el desempeño de las tuberías de ARN-seq en términos de predicción del desenlace de la enfermedad basada en la expresión génica. Estas muestras fueron proporcionadas por el Hospital de Niños de la Universidad de Colonia y secuenciadas en BGI utilizando la plataforma Illumina48. Las 176 muestras se tomaron de pacientes de riesgo alto que se definieron como aquellos con neuroblastoma en estadio 4 y edad > 18 meses o con tumores amplificados por MYCN de cualquier estadio o edad. El conjunto de datos SEQC-neuroblastoma se depositó en el Omnnibus de Expresión Génica con el número de adhesión GSE47792.

Predijimos dos variables clínicas: la supervivencia sin complicaciones (SSC), es decir, la aparición de acontecimientos como el progreso, la recaída o la muerte, y la supervivencia general (SG), es decir, la muerte. Para ambos criterios de valoración, los pacientes se dividieron en dos grupos (es decir, riesgos altos frente a riesgos bajos). Los pacientes de alto riesgo experimentaron un evento o murieron antes de un umbral de supervivencia predefinido, mientras que los pacientes de bajo riesgo experimentaron un evento o murieron después del umbral, o su último seguimiento superó el umbral. Los umbrales de supervivencia para la SSC y la SG fueron de dos y tres años, respectivamente. Los umbrales se eligieron para equilibrar el número de pacientes de alto y bajo riesgo. En el cuadro complementario S9 se proporcionan detalles del conjunto de datos SEQC-neuroblastoma.

También utilizamos un conjunto de datos ARN-seq de adenocarcinoma de pulmón de 87 muestras del repositorio del Atlas del Genoma del Cáncer (TCGA). El criterio de valoración de la predicción también fue la supervivencia, y utilizamos los mismos criterios para definir los grupos de alto y bajo riesgo con el umbral de supervivencia de dos años. El umbral de dos años se eligió para equilibrar el número de pacientes de alto y bajo riesgo. Los detalles del conjunto de datos TCGA-adenocarcinoma de pulmón se proporcionan en la Tabla suplementaria S10.

Filtrar el conjunto de datos de referencia de qPCR para producir un conjunto de genes de referencia

Debido a la variabilidad en las mediciones de qPCR y los desacuerdos entre las plataformas de qpcr7, filtramos el conjunto de datos de referencia de qPCR para retener genes que exhibían un comportamiento «correcto». A continuación, utilizamos estos genes para calcular las métricas de referencia (es decir, exactitud, precisión, fiabilidad y reproducibilidad). Tal proceso de filtrado se resume en la Fig. S1.

Comenzando con el conjunto inicial de 20,801 genes analizados con PrimePCR, filtramos estos genes para retener solo los genes que se cuantificaron como distintos de cero (es decir, detectados) y con valores de Ct (umbral de ciclo) ≤ 35 (35 indica la detección de una sola molécula en una muestra). El filtrado de los datos de PrimePCR resultó en 14.014 genes que también coincidieron con el transcriptoma AceView utilizado para mapear el conjunto de datos SEQC-benchmark RNA-seq.

Posteriormente, filtramos los 14,014 genes qPCR para retener solo 12,610 genes que exhibían el orden de titulación correcto (TO) y las relaciones de mezcla esperadas (EMR). Los detalles de este proceso se encuentran en la sección» Filtrado de genes qPCR por orden de titulación y proporciones de mezcla esperadas».

Por último, dado que algunas métricas de referencia, como la precisión y la precisión, son sensibles a genes de expresión cero o muy baja, seleccionamos genes que se expresaron como distintos de cero en todas las réplicas de todas las muestras de todos los sitios de secuenciación y las 278 tuberías de ARN – seq. El conjunto de referencia final contiene solo 10.222 genes qPCR (denominados «todos los genes») que se utilizaron para calcular las tres métricas de referencia para las tuberías de RNA-seq.

Según el estudio anterior, es más probable que los genes con menor expresión sean inconsistentes entre las tuberías 49. Por lo tanto, también identificamos un conjunto de genes de baja expresión en los 10,222 genes basados en la expresión promedio de qPCR de las muestras A, B, C y D. El 20% más bajo de los 10,222 genes (es decir, 2044 genes, denominados «genes de baja expresión») también se usaron para calcular el mismo conjunto de métricas de referencia para tuberías de ARN-seq. Este diseño nos permitió investigar la capacidad de las tuberías de ARN-seq para estimar la expresión génica de baja expresión.

Filtrado de genes qPCR por orden de titulación y relaciones de mezcla esperadas

Los conjuntos de datos de referencia SEQC (ARN-seq y qPCR) tienen propiedades únicas que permiten evaluar la corrección de la cuantificación. Después de identificar genes qPCR detectables (es decir, distintos de cero y Ct ≤ 35) y compatibles con AceView, utilizamos dos métricas (TO y EMR) para filtrar aún más el conjunto de datos qPCR de referencia, dejando solo genes qPCR «correctos». El EMR métricas de captura única mezcla de propiedades de los datos, es decir,

$$C= \frac{3}{4}A+\frac{1}{4}B\, \text{y }\,= \frac{1}{4}A+\frac{3}{4}B.$$

Debido a esta propiedad, todos los genes se espera que se expresa en una de las órdenes siguientes, dependiendo de la expresión relativa de las muestras a y B:

$$A\ge C\ge D\ge B \,\text{o }\, Un\le C\le D\le B.$$
$${\stackrel{-}{q}}_{s,\cdot ,k}=\frac{1}{N}\sum_{n=1}^{N}{q}_{s,n,k,}$$

el conjunto de qPCR de los genes que siga la correcta titulación orden

Para una sola replicar qPCR conjunto de datos (por ejemplo,, el conjunto de datos de PrimePCR que analizamos), la variabilidad inherente de una sola medición de qPCR puede resultar en algunos genes falsos negativos que siguen el correcto pero no se identifican. A partir de la literatura50,51, el coeficiente de variación para las mediciones replicadas de qPCR es generalmente del 15% o mayor, por lo que usamos este número para ajustar el margen para determinar si un gen sigue el A correcto. Matemáticamente, calculamos el rango de más y menos una desviación estándar de cada medición de qPCR y la usamos como margen. Las ecuaciones revisadas para \({K}_{A}\) son las siguientes:

$${K}_{A}={K}_{A,A\ge B}\cup {K}_{A,\le B}$$

donde \(a=1.15, b=0.85\)

Además, las muestras deben, además, presentan una específica relación de mezcla. Dado que la relación entre las muestras a y B es

$${R}_{A,B}=\frac{A}{B}$$

el EMR entre las muestras C y D,

$$EM{R}_{C,D}=\frac{3z\cdot {R}_{A,B}+1}{z\cdot {R}_{A,B}+3}\cdot \frac{z+3}{3z+1}$$
$${R}_{A,B}\in \left\equiv ,$$
$${R}_{C,D}\en \left\equiv \a la izquierda,\text{ y}$$
$$EM{R}_{C,D}\in \left\equiv ),$$

y, finalmente, determina un conjunto de genes que satisface la EMR criterio de la siguiente manera:

$${K}_{EMR}=\left\{k|\left({{R}_{C,D}^{Lower}\le {EMR}_{C,D}^{Superior}|}_{{k, R}_{C,D}\ge EM{R}_{C,D}}\right)\v \left({{R}_{C,D}^{Superior}\ge {EMR}_{C,D}^{Lower}|}_{{k, R}_{C,D}\le EM{R}_{C,D}}\right)\right\}$$

datos de RNA-seq tuberías de análisis de asignación, la cuantificación, y normalización

Hemos investigado 278 RNA-seq tuberías que incluye trece de la secuencia de asignación algorithms18,19,20,21,22,23,24,25,26,27,28,29, tres categorías de expresión de la cuantificación algorithms31,32,33, y siete expresión de métodos de normalización. Las tablas suplementarias S2-S4 resumen todas las opciones consideradas para cada componente de canalización (asignación de secuencias, cuantificación de expresiones y normalización de expresiones). Los trece algoritmos de mapeo investigados son Bowtie18, Bowtie219, BWA20, GSNAP21, Magic22 (una nueva tubería desarrollada por NCBI para el proyecto SEQC: ftp://ftp.ncbi.nlm.nih.gov/repository/acedb/Software/Magic), MapSplice23, Novoalign (un paquete comercializado desarrollado por Novocraft: https://www.novocraft.com/products/novoalign/), OSA24, RUM25, STAR26, Subread27, TopHat28 y WHAM29. Algunos usan mapeo sin empalmar de lecturas al transcriptoma, y otros realizan mapeo empalmado al genoma. Magic usa ambos en paralelo y compara la calidad de cada alineación para mantener la mejor en varios objetivos. Los algoritmos de asignación pueden informar solo de una asignación única o permitir varias ubicaciones de asignación por lectura. Los algoritmos de cuantificación incluyen métodos simples basados en recuentos (es decir, HTSeq31) y métodos probabilísticos basados en la distribución de Poisson aplicados a datos de mapeo genómicos (es decir, Cufflinks32) o transcriptómicos (es decir, RSEM33). La Magia, el RON y la Lectura secundaria (i. e., featureCounts52) las canalizaciones incluyen métodos de cuantificación integrados que entran en la categoría de métodos simples basados en recuento. Los métodos de normalización incluyen métodos de escalado simples (es decir , fragmentos por millón de fragmentos mapeados , fragmentos por kilobase de longitud de genes por millón de fragmentos mapeados, mediana y cuartil superior), métodos de escalado robustos (es decir, expresión de registro relativa y media recortada de valores m) y métodos integrados en tuberías específicas (es decir, índice de expresión mágica).

Asignación de secuencias

Asignamos secuencias a cada referencia en pasos sucesivos utilizando algoritmos de asignación sin empalmar o empalmados. La asignación sin empalmar se refiere a algoritmos que alinean secuencias de lectura completas (por ejemplo, Bowtie2, BWA y Novoalign), mientras que la asignación empalmada se refiere a algoritmos que dividen las lecturas en segmentos para acomodar espacios largos o intrones en una lectura (por ejemplo, TopHat y MapSplice). En el primer paso del mapeo sin empalmar, intentamos mapear todas las secuencias de extremos emparejados a la referencia ERCC/MT/rRNA (p. ej., Secuencias de Consorcio de controles de ARN externo, el genoma mitocondrial y secuencias de ARN ribosómico). Todos los pares de lectura no mapeados se asignaron al transcriptoma AceView. Finalmente, todos los pares leídos que no se asignaron a las referencias ERCC/MT/rRNA o AceView se asignaron a la referencia del genoma humano. Las coordenadas de mapeo transcriptómico se tradujeron a coordenadas de mapeo genómico y se fusionaron con los resultados de mapeo del genoma humano para producir los resultados finales (Fig. S21, panel izquierdo). Utilizamos Bowtie2 como mapeador para el primer paso de todas las tuberías de mapeo empalmadas(Fig. S21, panel derecho). Algoritmos de mapeo empalmados, ya sea que se asignen lecturas directamente al genoma humano (por ejemplo, MapSplice y GSNAP) o se asignen lecturas completas sin empalmar al transcriptoma y luego se fusionen estos resultados de mapeo con los resultados de mapeo empalmados de las lecturas restantes al genoma humano (por ejemplo, TopHat y OSA). La Tabla complementaria S2 resume todas las herramientas de mapeo investigadas en este estudio.

Bowtie2, GSNAP, Novoalign, TopHat y WHAM permiten controlar el número de asignaciones notificadas por par de lecturas. De forma predeterminada, estos algoritmos suelen informar de la mejor ubicación de asignación por par de lecturas. Sin embargo, algunos algoritmos de cuantificación pueden usar información sobre múltiples ubicaciones de mapeo ambiguas para mejorar la estimación de la expresión génica. Por lo tanto, además de los informes de una sola visita, también generamos resultados de mapeo que informaron hasta 200 visitas por lectura (múltiples visitas). También incluimos la tubería de mapeo de pajaritas con parámetros de mapeo específicos para la cuantificación con RSEM, como se describe en la siguiente sección33.

Las opciones de línea de comandos para todas las herramientas de alineación de secuencias se detallan en la Nota complementaria 1.

Cuantificación de la expresión génica

La etapa de cuantificación incluyó tres categorías de cuantificadores: cuantificadores basados en recuentos (es decir, HTSeq y cuantificadores incorporados para las tuberías Magic, RUM y Subread), cuantificadores basados en modelos de probabilidad para el mapeo genómico (es decir, , Gemelos) y cuantificadores basados en modelos de probabilidad para mapeo transcriptómico (es decir, RSEM). Las características clave de estos cuantificadores se resumen en el cuadro complementario S3. Cufflinks es un cuantificador basado en modelos de Poisson que estima las probabilidades de asignación de lectura en función de la información de alineación32. Es capaz de ensamblar transcripciones y cuantificar expresiones de genes o transcripciones. En este estudio, deshabilitamos la función de ensamblaje y proporcionamos el archivo GTF de anotación genómica como referencia de cuantificación. HTSeq es un cuantificador ingenuo basado en recuentos que asigna lecturas asignadas a genes31. HTSeq es capaz de cuantificar la expresión génica, pero no la expresión de transcripción. RSEM es también un cuantificador basado en modelos de Poisson que es similar en concepto a Cufflinks33. La información de las lecturas de múltiples visitas es importante tanto para los gemelos como para el RSEM. Estos algoritmos utilizan información de lectura de múltiples visitas para estimar con mayor precisión la expresión de genes o transcripciones.

Los resultados de mapeo de tuberías de alineación no siempre fueron compatibles con las tres categorías de cuantificadores. Los gemelos requieren que los resultados de alineación se ordenen por coordenadas de alineación y las lecturas de múltiples visitas se marquen con la etiqueta ‘ NH ‘ en el campo de atributos del archivo SAM. HTSeq requiere que los resultados de la alineación se ordenen por nombres leídos y que la etiqueta’ NH ‘ esté ausente en el archivo SAM. El RSEM solo cuantifica el mapeo transcriptómico, es decir, lee mapeados y reportados en coordenadas transcriptómicas. Además, RSEM solo maneja alineaciones sin separación. Por lo tanto, se requiere filtrado para eliminar las alineaciones con huecos. Debido a estos requisitos, procesamos previamente todos los resultados de alineación antes de la cuantificación. En resumen, veinte tuberías de alineación, incluidas las tuberías empalmadas, sin empalmar, de golpe único y de golpe múltiple, eran adecuadas para la cuantificación basada en recuento. Dieciséis tuberías de alineación eran adecuadas para gemelos, y solo diez eran adecuadas para RSEM. RSEM está diseñado específicamente para funcionar bien con pajarita. Por lo tanto, también incluimos esta canalización de mapeo y cuantificación incrustada.

Las opciones de línea de comandos para todos los instrumentos de cuantificación se detallan en la nota complementaria 1.

Normalización de la expresión génica

La normalización de los datos de ARN-seq permite la comparación entre muestras. Generalmente, los métodos de normalización corrigen el tamaño de la biblioteca (es decir, el número total de lecturas en una muestra), que es la fuente principal de varianza entre muestras. Investigamos siete métodos de normalización: fragmentos por millón de fragmentos mapeados (FPM), fragmentos por kilobase de longitud génica por millón de fragmentos mapeados (FPKM), mediana (Med.), cuartil superior (UQ), expresión logarítmica relativa (RLE), media recortada de valores M (TMM) e índice de expresión (EIndex, que es específico de la canalización Mágica) (véase la Tabla suplementaria S4). Describimos cada uno de estos métodos de normalización basados en la siguiente descripción matemática del conjunto de datos SEQC-benchmark.

$$\overline{x}_{s, \cdot ,k} = \frac{1}{N}\mathop \sum \limits_{n = 1}^{N} x_{s,n,k}$$

hemos definido el conjunto de genes presentes a ser

y el final de la presente conjunto de genes es

$$K_{p} = K_{p,BGI} \cap K_{p,PUEDE} .used

Utilizamos el mismo conjunto de gens presentes para todos los métodos de normalización para una tubería RNA-seq.

El número total de genes presentes para una muestra dada s y replicar n es

$$x_{s,n} = \mathop \sum \limits_{{k \in K_{p} }} x_{s,n,k} ,$$

y el promedio del número total de genes presentes para todos los datos de un único sitio está

$$\bar{x} = \frac{1}{4}\frac{1}{N}\mathop \sum \limits_{s} \mathop \sum \limits_{{n = 1}}^{N} x_{{s,n}}.$$

por Lo tanto, hemos definido FPM normalizada de expresión para cada muestra s, replicar n, y el gen k como

$$y_{s,n,k}^{DF} = \frac{{x_{s,n,k} \cdot \overline{x}}}{{x_{s,n} }}.$$

Mediana y alta cuartil normalizada de expresión para cada muestra s, replicar n, y el gen k se define como

$$y_{s,n,k}^{Med} = \frac{{x_{s,n,k} \cdot \tilde{x}}}{{\tilde{x}_{s,n} }}{\text{y }}y_{s,n,k}^{UQ} = \frac{{x_{s,n,k} \cdot \hat{x}}}{{\hat{x}_{s,n} }}.{}

Para la normalización FPKM, definimos la longitud de un gen k como \(\ell_{k}\), que es la longitud de la unión de todos los exones relacionados con el gen tal como se define en el transcriptoma AceView. La formulación original del FPKM utilizó arbitrariamente factores de escala de 1 × 103 para la longitud del gen y de 1 × 106 para el número total de fragmentos mapeados. Con el fin de mantener un rango dinámico comparable entre todos los métodos de normalización, en su lugar, escalamos por la longitud media del gen y el recuento total promedio de todos los genes presentes. La longitud media de todos los presentes de los genes

$$\overline{\ell } = \frac{1}{{\left| {K_{p} } \right|}}\mathop \sum \limits_{{k \in K_{p} }} \ell_{k} .$$

por Lo tanto, reescalado FPKM normalizada de expresión para cada muestra s, replicar n, y el gen k es

$$y_{s,n,k}^{FPKM} = \frac{{x_{s,n,k} \cdot \overline{\ell } \cdot \overline{x}}}{{x_{s,n} \cdot \ell_{k} }}.

Los métodos de normalización TMM y RLE son similares a la normalización FPM, pero introducen un factor de escala adicional para ajustar el tamaño de la biblioteca. Utilizamos el paquete bordeador en R para estimar un factor de escala para cada replicante de muestra36,53. El método TMM selecciona una biblioteca de referencia de un conjunto de bibliotecas replicadas de muestra y, a continuación, calcula las relaciones de expresión logarítmica gene-wise (valores M) y los valores medios de expresión logarítmica gene-wise (valores A) entre la biblioteca de destino y la biblioteca de referencia. Los números extremos en los valores M y A se recortan, y el factor de escala para la biblioteca de destino es el promedio ponderado de los valores M restantes. El método RLE determina un factor de escala definiendo primero la biblioteca mediana como la media geométrica gene-wise de los replicados de la muestra35. La relación mediana de cada biblioteca de destino con la biblioteca mediana se toma como factor de escala. Las expresiones normalizadas de TMM y RLE para cada muestra s, replicar n y gen k se definen como:

donde \(\hat{f}_{s, n}^{TMM}\) y \(\hat{f}_{s,n}^{RLE}\) son el factor de escala para la muestra s,replicar n.

Métricas de rendimiento de tuberías de RNA – seq

Las métricas de referencia para tuberías de RNA-seq se resumen en la Tabla suplementaria S7.

Precisión medida como desviación de las referencias qPCR

$ $ {\stackrel{-}{y}}_{s,\cdot ,k}=\frac{1}{N}\sum_{n=1}^{N}{y}_{s,n,k} Given

Dadas las muestras A y B, la desviación absoluta de la relación logarítmica de la expresión basada en ARN-seq de la expresión basada en qPCR la expresión para un gen k es

$ $ \Delta_{\frac{A}{B},k} = \left | \log_2\left ( \frac{\bar{x}_{A,., k}} {\bar {x}_{B,.,k}} \ right) – \log_2 \ left (\frac {\bar {q} _ {A,., k}} {\bar {q} _ {B,.,k}} \ right) \right|,

y la métrica de precisión final se definió como la mediana de todos los \({\Delta} _{{\frac{A}{B},k}}\), \(k = 1\ldots K\).

Precisión medida como variación de la expresión génica en bibliotecas replicadas

Calculamos el coeficiente de variación (CoV) para cada gen y cada muestra en cuatro bibliotecas replicadas de la siguiente manera:

$$CoV_{s,k} = \frac{{sd\left( {x_{s, \cdot ,k} } \right)}}{{\overline{x}_{s, \cdot ,k} }},$$

la Fiabilidad se mide como intra-muestra la correlación de la expresión de genes

La fiabilidad de un sistema de medición puede ser evaluada mediante el coeficiente de correlación intraclase (ICC)54,55. ICC es aplicable a mediciones que se pueden organizar en grupos, y describe cómo las mediciones del mismo grupo son similares entre sí. La definición ICC moderna toma prestado el marco de análisis de varianza (ANOVA), o más específicamente ANOVA con efectos aleatorios55. El tipo de ANOVA depende del diseño experimental y generalmente sigue la definición del artículo de Shrout publicado en 197955. ICC(1,1) y ICC (1,k) se basan en el modelo de efectos aleatorios unidireccionales y son aplicables en el caso de que cada grupo sea evaluado por un conjunto diferente de evaluadores k seleccionados aleatoriamente de una población mayor de evaluadores. ICC(2,1) y ICC (2,k) se basan en el modelo de efectos aleatorios bidireccionales y son aplicables al caso de que se preseleccione una muestra aleatoria de evaluadores k de una población mayor y cada evaluador evalúe a cada grupo exactamente una vez (p. ej., cada evaluador evalúa n grupos en conjunto). CIC (3,1) y CIC(3,k) se basan en el modelo de efectos mixtos bidireccionales y son aplicables al caso de que cada grupo sea evaluado por cada uno de los mismos evaluadores k, que son los únicos evaluadores de la población. El segundo parámetro en ICC(,) indica si el ICC debe medir la fiabilidad de una sola medición o el promedio de las mediciones k.

Para el conjunto de datos de referencia SEQC con bibliotecas replicadas para cada muestra, ICC(1,1) o ICC(1,k) se ajustaron a nuestro objetivo, ya que, para un gen g, la expresión génica de bibliotecas replicadas para diferentes muestras (o grupos diferentes en el contexto anterior) no se evaluó exactamente en las mismas condiciones (o evaluada por los mismos evaluadores en el contexto anterior). Elegimos usar ICC (1, k) ya que las bibliotecas replicadas están disponibles para la mayoría de los experimentos. Matemáticamente, un modelo de efectos aleatorios unidireccionales se puede formular como

$ $ Y_{ij} = \mu + \alpha_{j} +{\varepsilon_{ij}},
di ICC\left( {1,k} \right) = \frac{BMS-WMS}{{BMS}}, calculated

Calculamos ICC para cada gen k, \(k = 1 \ldots K\), y luego se utilizó la mediana de todos los CI como medida final de fiabilidad.

También hemos investigado otras métricas potenciales, como la reproducibilidad, que se define como la correlación de Spearman entre dos bibliotecas replicadas de la misma muestra (Nota Complementaria 2). La correlación de Spearman varió de 0.993 a 0,996 (Fig. S8) usando AllGenes. Descartamos la métrica de reproducibilidad debido al rango dinámico relativamente pequeño.

Evaluación de la utilidad de las métricas de referencia para la selección de tuberías de RNA-Seq

Clasificamos las tuberías de RNA-seq en función del rango promedio de las tres métricas de referencia (es decir, precisión, precisión y fiabilidad). A continuación, evaluamos la utilidad de las métricas de referencia examinando si las canalizaciones de buen y mal rendimiento identificadas con base en las métricas de referencia eran informativas para inferir el rendimiento de la predicción del desenlace de la enfermedad basada en la expresión génica y la significación estadística de la estratificación de los pacientes para todas las variables clínicas (es decir, las variables SEQC-SSC y SG del neuroblastoma y la variable TCGA-supervivencia del adenocarcinoma de pulmón).

En primer lugar, para las 278 tuberías representativas de ARN-seq aplicadas al conjunto de datos de referencia SEQC, calculamos el rango promedio utilizando un subconjunto de las métricas de referencia como indicador de rendimiento final para cada tubería. En total, teníamos 6 métricas (3 métricas de referencia × 2 conjuntos de genes) e investigamos 12 subconjuntos (4 × 3) de las 6 métricas utilizando los siguientes criterios:

  1. (1)

    Cuatro combinaciones de las tres métricas de referencia con al menos dos en un subconjunto: una combinación con las tres métricas de referencia, tres combinaciones con dos de las tres métricas de referencia.

  2. (2)

    Tres subconjuntos formados por métricas derivadas de todos los genes, las derivadas de genes de baja expresión, o una combinación de ambos.

En segundo lugar, para cada una de las 278 tuberías representativas de ARN-seq (156 para el criterio de valoración de supervivencia de TCGA-adenocarcinoma de pulmón), calculamos valores de AUC y CCM anidados de validación cruzada, como se describe en la sección «Método» «Modelos predictivos de neuroblastoma y adenocarcinoma de pulmón», lo que resultó en 834 (468 para el criterio de valoración de supervivencia de TCGA-adenocarcinoma de pulmón) valores de AUC y CCM para cada criterio de valoración clínico (p. ej., 278 clasificadores de tuberías × 3, o 156 clasificadores de tuberías × 3) (Cuadros complementarios S11,S12). También modelamos las funciones de supervivencia utilizando el análisis de Kaplan–Meier para cada canalización, como se describe en la sección «Método» «Análisis de supervivencia de Kaplan–Meier». Para cada canalización de ARN-seq, resumimos el rendimiento de la predicción del desenlace de la enfermedad basada en la expresión génica utilizando tanto el AUC promedio como el CCM entre los clasificadores y la tasa de éxito de la estratificación de los pacientes (p. ej., separación estadísticamente significativa de dos curvas de Kaplan–Meier) en todas las iteraciones y clasificadores en el marco de validación cruzada anidado.

Finalmente, identificamos el 10% superior de tuberías de buen rendimiento y el 10% inferior de tuberías de mal rendimiento en función del rango promedio de un subconjunto de las tres métricas de referencia. El rendimiento de predicción correspondiente (p. ej., AUC y MCC) de los oleoductos de buen rendimiento se compararon con los oleoductos de mal rendimiento utilizando la prueba de suma de rango de Wilcoxon unilateral con la hipótesis nula de que la mediana del primer grupo no era mayor que la del segundo grupo.

Modelos predictivos de neuroblastoma y adenocarcinoma de pulmón

Evaluamos el rendimiento de 278 tuberías de ARN-seq en términos de toma de decisiones basada en la expresión génica utilizando el conjunto de datos de neuroblastoma SEQC48. El conjunto de datos SEQC-neuroblastoma y los criterios de valoración clínicos relacionados se resumen en la Tabla complementaria S9. Las tuberías de ARN-seq se evaluaron en términos de predicción de los desenlaces de pacientes con neuroblastoma para dos criterios de valoración clínicos mediante validación cruzada anidada (Fig.Suplementaria. S13) 56,57. También evaluamos de manera similar el rendimiento de 156 tuberías de ARN-seq aplicadas al conjunto de datos TCGA-adenocarcinoma de pulmón para predecir el desenlace de la enfermedad. El conjunto de datos TCGA-adenocarcinoma de pulmón y los criterios de valoración clínicos asociados se resumen en la Tabla complementaria S10.

La validación cruzada anidada implica entrenamiento y pruebas de un modelo de predicción óptimo. Esto se logra utilizando la validación cruzada interna o de optimización triple, aplicada al subconjunto de entrenamiento de la validación cruzada externa quíntuple. Una vez identificados los parámetros finales del modelo de predicción óptima (es decir, los hiperparámetros clasificadores y el tamaño de la entidad), el modelo final se entrena utilizando todo el subconjunto de entrenamiento y, a continuación, se prueba utilizando el pliegue restante de la validación cruzada exterior quíntuple. Este proceso se repitió durante diez iteraciones. Realizamos una validación cruzada anidada por separado para cada uno de los tres clasificadores (p. ej., refuerzo adaptativo, regresión logística y máquinas de vectores de soporte) y utilizó el método de selección de entidades de redundancia mínima y relevancia máxima (mRMR) para elegir tamaños de entidades óptimos dentro del rango de 5 a 40 con el tamaño de paso de 558.

Análisis de supervivencia de Kaplan–Meier

Para cada canalización y clasificador de ARN-seq (es decir, 278 clasificadores de canalizaciones × 3 para las variables SEQC-neuroblastoma y 156 clasificadores de canalizaciones × 3 para la variable TCGA-adenocarcinoma de pulmón), modelamos las funciones de supervivencia de Kaplan-Meier en función de las etiquetas previstas de cada muestra. A continuación, se utilizó la prueba de rango logarítmico de dos colas para determinar si las curvas de supervivencia estimadas para cada grupo de pacientes previsto eran estadísticamente diferentes.

Análisis de varianza y cálculo de la contribución de cada factor de canalización ARN-seq a la varianza general de la canalización

Utilizamos el análisis de varianza (ANOVA) para determinar si cada factor de canalización ARN-seq contribuye significativamente a la varianza de cada una de las tres métricas de referencia (es decir, exactitud, precisión y fiabilidad), así como a la varianza del rendimiento de predicción (es decir, AUC y MCC). Para cada una de las tres métricas de referencia, utilizamos un modelo lineal (función R «lm») para ajustar los datos de las 278 tuberías utilizando la métrica como variable dependiente y los factores de tubería RNA-seq como variables categóricas independientes. Consideramos los siguientes factores como variables categóricas independientes: algoritmo de mapeo, estrategia de mapeo (es decir, empalmado vs.no empalmado), informes de mapeo (es decir, golpe único vs. golpe múltiple), algoritmo de cuantificación y algoritmo de normalización. Se incluyeron todos los factores y sus interacciones bidireccionales en el modelo lineal. Para cada uno de los puntos finales de predicción, aplicamos la misma técnica para ajustar los datos de las 278 tuberías utilizando el AUC promedio o el CCM como variable dependiente y el mismo conjunto de factores de tubería ARN-seq como variables categóricas independientes. Luego realizamos el ANOVA en el modelo lineal (función R «anova»). ANOVA calcula una «suma de cuadrados» (es decir, varianza) atribuida a cada factor o interacción y utiliza una prueba F para determinar si la varianza es estadísticamente significativa. Calculamos el porcentaje que cada factor o interacción contribuye a la varianza total calculando la relación de «suma de cuadrados» para cada factor a la suma total de cuadrados.

Análisis de regresión

Se investigó la relación entre los perfiles de alineación o las características de distribución de expresión génica y las métricas de referencia. Los perfiles de alineación incluyeron el número total de fragmentos mapeados, el número total de lecturas que abarcan la región intrónica, el número total de lecturas con inserciones o eliminaciones, el número total de lecturas perfectamente coincidentes, el número total de lecturas con como máximo un desajuste y el número de desajustes por lectura mapeada. Cada algoritmo de alineación estaba representado por las estadísticas promedio de 2 sitios de secuenciación, 4 muestras, 4 bibliotecas replicadas y 2 carriles. Utilizando el paquete de «MASA» en R, adoptamos el método de estimación M con ponderación Huber para ajustar modelos de regresión lineal robustos entre una variable dependiente (rendimiento de la métrica de referencia) y una variable explicativa (un perfil de alineación). El método de estimación M con ponderación Huber es un método de regresión robusto en presencia de valores atípicos. Las características de distribución de expresión génica incluyeron el cuartil inferior, mediana, cuartil superior, máximo, rango intercuartílico, desviación estándar, asimetría, curtosis y entropía de una distribución de expresión génica. Utilizamos el mismo método de estimación M con ponderación Huber para ajustar modelos de regresión lineal robustos, y luego reportamos el error estándar residual para cada modelo.

Descargo de responsabilidad

Las opiniones presentadas en este artículo no reflejan necesariamente la opinión o política actual o futura de la Administración de Alimentos y Medicamentos de los Estados Unidos. Cualquier mención de productos comerciales es para aclaración y no pretende ser un respaldo.

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *