Estimación del mapeo de frecuencia y asociación de alelos utilizando datos de secuenciación de próxima generación

El alelo menor es el alelo menos frecuente en la población de un sitio variable. En primer lugar, describimos dos enfoques principales para estimar la frecuencia de alelos menores (MAF) en un sitio particular del genoma. El primer enfoque implica inferir genotipos individuales y tratar los genotipos inferidos como completamente precisos al estimar el MAF. A continuación, examinamos el rendimiento de un marco de verosimilitud que tiene en cuenta directamente la incertidumbre en la asignación de genotipos. A lo largo de nuestro trabajo, asumimos que todos los sitios segregados son bialélicos.

Estimación de MAF a partir de genotipos llamados

Una forma de estimar el MAF a partir de datos de secuenciación de próxima generación es llamar primero a un genotipo para cada individuo utilizando datos de secuenciación, y luego usar esos genotipos como si fueran los verdaderos. Este fue el enfoque utilizado tradicionalmente para los datos de genotipos y los datos de secuenciación de Sanger. No está claro qué tan bien funcionará cuando se aplique a datos de secuenciación de próxima generación.

Se puede utilizar un enfoque de máxima verosimilitud para inferir el genotipo de cada individuo a partir de los datos de secuenciación de próxima generación. En cada sitio j, para cada individuo i, la probabilidad para cada uno de los tres genotipos posibles (suponiendo que conocemos el alelo menor) se da como:

(1)

donde D i,j son los datos de secuenciación observados en el individuo i en el sitio j, g i , j ∈ {0, 1, 2} es el número de alelos menores contenidos en el genotipo de cada individuo, yy control para errores de secuenciación y calidades base de lectura, respectivamente. Los datos de secuenciación observados para cada individuo se pueden considerar como la alineación de las lecturas en el sitio j teniendo en cuenta las puntuaciones de calidad de lectura. Esto se representa como la probabilidad de genotipo y se encuentra en el archivo de probabilidad de genotipo (GLF) que se produce en muchos programas que analizan datos de secuenciación de próxima generación, como SOAPsnp y MAQ .

Para asignar un genotipo a un individuo en particular, se puede calcular la probabilidad de cada uno de los tres genotipos posibles para el individuo. Entonces se puede asignar el genotipo con mayor probabilidad. Sin embargo, los investigadores a menudo prefieren un criterio de llamada más estricto y no asignarán un genotipo a un individuo a menos que el genotipo más probable sea sustancialmente más probable que el segundo más probable. Aquí los tres genotipos posibles se ordenan por sus probabilidades: , donde g(k)corresponde al genotipo con la k-ésima mayor probabilidad. Con un umbral dado f, se puede llamar al genotipo g (1) if . De lo contrario, no se llama a un genotipo y se considera que falta el genotipo del individuo. Un valor umbral común de f es 1, lo que indica que el genotipo más probable es al menos 10 veces más probable que el segundo genotipo más probable. Tenga en cuenta que este tipo de filtrado puede resultar en una mayor confianza para el genotipo «llamado», pero también resulta en más datos faltantes.

Estimador de máxima verosimilitud de la frecuencia de alelos

En lugar de estimar el MAF a partir de los llamados genotipos, un método de máxima verosimilitud (ML) introducido por Kim et al. (véase también Lynch para un enfoque similar) estima directamente los FMA y tiene en cuenta la incertidumbre del genotipo. Específicamente, dado un alelo menor, la probabilidad de observar los datos de secuencia en cada individuo i se obtiene sumando las probabilidades correspondientes a los tres genotipos posibles.

Supongamos que las tres probabilidades de genotipo definidas en la Ecuación 1 están disponibles. Usando la misma notación que la anterior, sean D j y p j los datos de secuenciación observados en el sitio j y el MAF correspondiente, respectivamente. La probabilidad de genotipo dada la frecuencia de alelos menores se puede calcular asumiendo el equilibrio de Hardy-Weinberg (HWE). Entonces, asumiendo la independencia entre individuos, la probabilidad del MAF en este locus es un producto de todas las probabilidades calculadas a través de todos los N individuos:

(2)

La estimación de ML de p j se puede calcular maximizando directamente la probabilidad de un espacio de parámetros restringido Método de Shanno (BFGS) o mediante el algoritmo de maximización de expectativas (EM). Cuando se utiliza el algoritmo EM, se calcula la expectativa posterior de un genotipo para cada individuo, y la media de esos posteriores se actualiza repetidamente. Nuestra implementación de BFGS fue más rápida que el algoritmo EM. Por ejemplo, para obtener estimaciones de 100,000 sitios, los BFG tomaron ~16 segundos, pero los EM tomaron ~100 segundos. Sin embargo, la diferencia en la velocidad puede ser específica de la implementación. En nuestro caso, para ambos métodos, dejamos de actualizar los parámetros cuando el aumento de la probabilidad era inferior a 0,001.

Estimador de máxima verosimilitud con alelo menor incierto

En la práctica, a menudo el segundo nucleótido más común entre los individuos se puede usar como alelo menor. Sin embargo, para SNPs raros (por ejemplo, MAF < 1%), es difícil determinar qué alelo es el alelo menor, ya que los cuatro nucleótidos pueden aparecer en algunas lecturas debido a errores de secuenciación. Para hacer frente a esta situación, ahora describimos un marco de probabilidad que tiene en cuenta la incertidumbre en la determinación del alelo menor.

Supongamos que para el sitio j conocemos el alelo mayor M. Tenga en cuenta que decidir cuál de los dos alelos comunes es probable que sea el mayor no es importante, ya que nos preocupa principalmente estimar las frecuencias en SNP raros. Además, para alelos con frecuencias intermedias (alrededor del 50%), la distinción entre alelos mayor y menor es menos importante. Asigne los otros tres nucleótidos no principales m1, m2 y m3. La verosimilitud introducida en la Ecuación 2 asume un alelo mayor fijo M y un alelo menor fijo m. Por lo tanto, para permitir la incertidumbre en la designación del alelo menor, la función de verosimilitud se puede modificar como:

(3)

Además, suponiendo que cualquiera de los tres posibles menor alelos es igualmente probable, obtenemos:

(4)

donde . Dado que puede ser muy pequeño con grandes conjuntos de datos (por ejemplo, con muchos individuos), es útil calcular la probabilidad en la escala logarítmica. Ordene las tres probabilidades de registro condicionales como (l(1), l(2), l(3)), donde l(1) es la mayor. Luego,

Prueba G usando genotipos llamados para mapeo de asociación

En estudios de asociación, se dice que los SNP que muestran diferencias significativas en la frecuencia de los alelos entre los casos y los controles están asociados con el fenotipo de interés. El mapeo de asociaciones se puede realizar utilizando datos de estudios de secuenciación de próxima generación. Primero discutimos los enfoques que requieren llamar a genotipos individuales y luego realizamos una prueba de asociación utilizando los genotipos llamados. En este enfoque, se llama primero a un genotipo para cada individuo. Los genotipos pueden filtrarse o no. Asumiendo la independencia entre individuos y HWE, se puede construir una tabla de contingencia de 2 × 2 contando el número de alelos mayores y menores tanto en los casos como en los controles. Esto conduce a la conocida prueba de razón de verosimilitud para la independencia, la prueba G:

(5)

donde O k,h es la frecuencia observada en una celda, y E k,h es la frecuencia esperada bajo la hipótesis nula en la que la frecuencia alélica es la misma entre casos y controles. La conocida prueba chi-cuadrado de Pearson es asintóticamente equivalente a la prueba G. Si la tabla se genera a partir de genotipos verdaderos, el estadístico G sigue asintóticamente una distribución chi-cuadrado con 1 grado de libertad (χ2(1)). Sin embargo, en nuestros estudios, construimos el estadístico G utilizando genotipos «llamados», por lo que es posible que el HWE no se mantenga debido a la sobre-y sub – llamada de los heterocigotos. Además, la construcción de la estadística de la prueba contando genotipos» llamados » en lugar de genotipos «observados» probablemente introduce variabilidad adicional. Por lo tanto, la teoría estadística puede no ser válida más. Tenga en cuenta que cuando no se llama a un genotipo para un individuo determinado, los datos se consideran faltantes y no se incluyen en la tabla de 2 × 2.

Prueba de razón de verosimilitud que tiene en cuenta la incertidumbre en los genotipos observados para el mapeo de asociación

En lugar de llamar a los genotipos, el marco de verosimilitud permite la incertidumbre en los genotipos y las pruebas en cada sitio j si la frecuencia de los alelos es la misma entre los casos y los controles. Formalmente, calculamos la probabilidad de las hipótesis H O: p j, 1 = p j, 2 (=p j ,0) y H A : p j ,1 ≠ p j ,2, donde p j ,1 y p j ,2 son de la MAFs en casos y controles, respectivamente.

Suponiendo que se conocen los alelos menores (m) y mayores (M), la probabilidad de la frecuencia del alelo menor se puede calcular como se describe en la Ecuación 2, y la estadística de la prueba de razón de verosimilitud se calcula como:

(6)

donde y son los datos observados para los casos y controles, respectivamente, y y son los Emv de la MAFs en casos y controles, respectivamente.

Si el alelo menor es desconocido, la probabilidad bajo la hipótesis nula se calcula como en la Ecuación 3, y la estadística de TRL se modifica como:

(7)

donde D j es los datos observados para ambos casos y controles, y es la frecuencia del alelo bajo la hipótesis nula. Otras notaciones son las mismas que en la ecuación 6.

Estimación de MAF en datos simulados

Comparamos las estimaciones de la frecuencia de los alelos en datos simulados utilizando genotipos verdaderos (True), llamados genotipos sin filtrado (Llamada NF), llamados genotipos con filtrado (f = 1; Llamada F) y el método de máxima verosimilitud (ML). Para SNPs raros, el tipo de alelo menor a menudo no es aparente. Al llamar a genotipos, se supone que el segundo nucleótido más común es el alelo menor. El método ML incorpora directamente la incertidumbre en la determinación del alelo menor y, a menos que se indique lo contrario, se muestran los resultados utilizando el método del alelo menor desconocido (Ecuación 3). Tenga en cuenta que el método alelo menor desconocido ML funciona de manera similar al método alelo menor conocido ML, pero el primero es mejor para SNPs muy raros (archivo adicional 1).

Primero evaluamos qué tan bien los diferentes enfoques fueron capaces de estimar el MAF en 200 individuos en un rango de profundidades de secuenciación para 1,000 SNPs con un MAF verdadero del 5%. La Figura 1 muestra cuadros de distribución de los FMA estimados utilizando los cuatro enfoques diferentes. Como se esperaba, para datos de mayor cobertura, como una profundidad individual de 12×, todos los métodos funcionan tan bien como cuando los genotipos se conocen con certeza (Verdadero). Sin embargo, cuando la profundidad disminuye, las estimaciones del MAF obtenidas por los genotipos de primera llamada se sesgan. Por ejemplo, la mediana de MAF estimada utilizando el método Call F es de 5,3% a 6× cobertura y de 12,5% a 2×. La razón del sesgo ascendente es que se hace más difícil llamar heterocigotos, ya que los heterocigotos verdaderos a menudo parecen errores de secuenciación. Por lo tanto, más heterocigotos que homocigotos menores tienden a tener genotipos faltantes. Sin embargo, el sesgo general en las estimaciones de MAF de los llamados genotipos no siempre está en una dirección (no se muestran los datos). Curiosamente, el sesgo parece ser peor para el método de llamada F que para el método de llamada NF. Este patrón puede parecer contrario a la intuición, ya que filtrar las llamadas del genotipo parece disminuir la probabilidad de llamar heterocigoto a un error de secuenciación. Sin embargo, el método de llamada F también resulta en una mayor cantidad de datos faltantes, ya que muchos homocigotos para el alelo mayor no se llamarán debido a errores de secuenciación. Por lo tanto, en este caso, llamar genotipos sin filtrar parece ser la mejor estrategia que filtrar genotipos cuando se trata de estimar el MAF.

Figura 1
figura 1

Estimaciones de la frecuencia de los alelos en sitios con un MAF verdadero del 5% para diferentes profundidades de cobertura. En cada profundidad, se simularon 1000 sitios utilizando 200 individuos, y en cada sitio, se calcula una estimación de la frecuencia de los alelos utilizando: (1) genotipos verdaderos (True); (2) llamados genotipos sin filtrado (Llamada NF); (3) llamados genotipos con filtrado (Llamada F); y (4) el método de máxima verosimilitud (ML). Para obtener más información sobre los métodos de estimación, consulte Métodos.

Los resultados son radicalmente diferentes para el nuevo ML método. Este método proporciona estimaciones imparciales del MAF (mediana de ~4,9%) en un rango de profundidades. Incluso a 2×, las estimaciones muestran solo una varianza ligeramente mayor que las basadas en los genotipos verdaderos.

También comparamos el error cuadrado medio estimado (MSE; Expectativa () de las diferentes estimaciones del MAF en un rango de profundidades de secuenciación (Figura 2). El método ML tiene un MSE más bajo que los métodos de llamada con 50 o 200 individuos. En particular, el MSE calculado basado en el método de Llamada F es mucho más alto que los de los otros métodos, especialmente cuando la profundidad disminuye. El MSE de las estimaciones del MAF basadas en los genotipos verdaderos refleja el límite inferior del MSE y no es constante a través de profundidades debido a la varianza de muestreo y a un tamaño de muestra finito. Utilizando 50 individuos, el MSE se acerca a 0,0005 con una profundidad creciente y cuando se utiliza un tamaño de muestra de 200 individuos, se acerca a 0,0013 con una profundidad creciente.

Figura 2
figura 2

Error squred medio (MSE; Esperado ) de cuatro tipos de estimadores de frecuencia de alelos para diferentes tamaños de muestra (panel izquierdo y derecho) y profundidades de cobertura (eje x). En cada profundidad, el MSE se calculó a partir de las estimaciones de la frecuencia de los alelos realizadas utilizando cuatro métodos diferentes: True, Call NF, Call F y ML (para más detalles de los métodos, consulte la leyenda de la Figura 1).

en General, el nuevo ML método realiza fuera de genotipo llamar a los métodos.

Estimación de una distribución de MAF a partir de datos simulados

A continuación, examinamos cómo se realizaron los diferentes enfoques de estimación para estimar la proporción de SNP a diferentes frecuencias en la población (similar al espectro de frecuencias del sitio, pero basado en la frecuencia del alelo de la población en lugar de la frecuencia de la muestra). Aquí simulamos 20,000 SNPs donde la distribución de los MAF verdaderos siguió la distribución estacionaria estándar para un tamaño de población efectivo de 10,000 (ver Métodos). Tenga en cuenta que en la práctica, sin embargo, es muy difícil distinguir un SNP muy raro de un error de secuenciación. Por lo tanto, para fines de comparación con datos reales, descartamos SNP con MAF estimado inferior al 2%. La Figura 3 muestra la proporción de SNP que caen en cada compartimiento de frecuencia diferente después de excluir aquellos SNP con MAF estimado< 2%.

Figura 3
figura 3

Distribución de frecuencias alélicas de SNPs simuladas asumiendo la distribución estacionaria estándar de frecuencias alélicas. En cada profundidad (cada panel), se simularon 20.000 SNP, y para cada SNP, se obtuvieron estimaciones del MAF utilizando cuatro métodos diferentes (ver el pie de foto de la Figura 1). Luego, para cada método (cada color), solo se utilizan sitios con frecuencias alélicas estimadas > 2% para generar cada histograma (eje x).

Como se esperaba, con una alta profundidad de cobertura, como 10× por individuo, todos los métodos proporcionan distribuciones de MAF estimadas que son similares a la distribución esperada basada en los genotipos verdaderos (Figura 3). Con una profundidad de cobertura menor, como menos de 4× por individuo, las distribuciones de MAF obtenidas por métodos de llamada de genotipos se apartan significativamente de la distribución de MAF esperada basada en genotipos verdaderos (Figura 3). En particular, estos métodos sobreestiman la proporción de SNP de baja frecuencia. Por ejemplo, la proporción esperada de SNP en el segundo bin (MAF estimado entre 2-4%) es del 18%. La proporción correspondiente basada en el método de llamada NF a una profundidad de 4× es del 26%, que es 1,4 veces mayor de lo esperado. La sobreestimación de la proporción de SNP de baja frecuencia se produce debido a la confusión de errores de secuenciación con heterocigotos verdaderos, lo que resulta en genotipos heterocigotos en sobrecalificación. La magnitud de esta inflación difiere según los diferentes cortes de filtrado, pero un corte más grande no necesariamente aumenta o disminuye la inflación.

La imagen es completamente diferente para el método ML. La distribución estimada de MAF obtenida del nuevo método de ML sigue de cerca la distribución verdadera incluso con profundidades de cobertura poco profundas. Aquí casi no hay exceso de SNPs de baja frecuencia. A una profundidad de 4×, la proporción de SNPs en la segunda bandeja del histograma es del 18,4%, que está muy cerca de la proporción esperada (18%). Por lo tanto, se pueden hacer estimaciones más confiables del espectro de frecuencias a partir de datos de baja cobertura utilizando nuestro enfoque de probabilidad que utilizando los enfoques de llamada de genotipo.

Mapeo de asociación en datos simulados

Comparamos el rendimiento de los métodos que tratan genotipos inferidos como genotipos verdaderos en pruebas de asociación (utilizando una prueba G) con nuestra prueba de razón de verosimilitud (TRL) que tiene en cuenta la incertidumbre en los genotipos. Examinamos la distribución de la prueba estadística bajo la hipótesis nula de diferencia de frecuencia sin alelos entre casos y controles. También comparamos el poder de los diferentes enfoques.

Con tamaños de muestra razonablemente grandes, la teoría asintótica estándar sugiere que bajo la hipótesis nula, tanto el estadístico G como el estadístico LRT siguen una distribución chi-cuadrado con un grado de libertad (χ2(1)). Por lo tanto, hemos comparado la distribución nula de la estadística G calculada sobre la base de métodos de llamada, así como la estadística LRT, con la distribución χ2(1) utilizando gráficas QQ (Figura 4). Simulamos 5000 SNPs a través de una variedad de profundidades de secuenciación en 500 casos y controles, donde el MAF utilizado para simular genotipos fue del 5% en ambos casos y controles. La distribución del estadístico G calculada utilizando los genotipos verdaderos muestra una muy buena correspondencia con una distribución χ2(1). Sin embargo, la distribución del estadístico G calculada en base a los llamados genotipos se aparta sustancialmente de una distribución χ2(1). Llamar a los genotipos y luego tratarlos como precisos produce un gran exceso de señales falsas positivas si los valores de p se calculan utilizando una distribución χ2(1). Por ejemplo, a una profundidad de 2×, el 11% de los SNP tenían un valor de p inferior al 5%, en comparación con el 5% esperado. El efecto es causado por un aumento en la varianza, debido al sobrecalentamiento de homocigotos como heterocigotos, en la prueba alélica utilizada aquí para detectar asociación. Las pruebas genotípicas como la prueba de tendencia de Armitage, que son robustas a las desviaciones del equilibrio de Hardy-Weinberg, no muestran un aumento similar en la tasa de falsos positivos (archivo adicional 2). Consistente con esta observación, el filtrado de los llamados genotipos resulta en una disminución en la fracción de pruebas significativas cuando se usa la prueba G, aunque el filtrado no resuelve completamente el problema. Por otro lado, la estadística de TRL muestra solo una desviación muy leve de una distribución χ2(1) para profundidades de cobertura de 2× o 5×.

Figura 4
figura 4

QQ-gráficas que comparan la distribución nula del estadístico de prueba de interés con una distribución χ2(1). Cada columna corresponde a una estadística de prueba diferente: (1) Estadística G calculada utilizando los genotipos verdaderos (True); (2)estadística G calculada utilizando genotipos llamados sin filtrado (Llamada NF); (3) Estadística G calculada utilizando genotipos llamados con filtrado (Llamada F); y (4) estadística de prueba de razón de verosimilitud con alelo menor desconocido (TRL). Suponiendo 500 casos y 500 controles, bajo la hipótesis nula, se simuló un conjunto de 5.000 sitios con un MAF del 5% con una profundidad de secuenciación de 2× (paneles superiores) y 5× (paneles inferiores). El factor «Inflación» se muestra en la esquina superior izquierda de cada figura.

también genera receiver operating characteristic (ROC) para cada una de las diferentes pruebas de asociación. Estas curvas muestran la potencia de la prueba a diferentes tasas de falsos positivos. Dado que las distribuciones de algunas de las estadísticas de prueba no siguen la distribución χ2(1) bajo la hipótesis nula, para hacer una comparación justa, obtuvimos el valor crítico para cada tasa de falsos positivos basado en la distribución nula empírica. La potencia se calcula como la fracción de los loci de enfermedad simulados que tienen una estadística que excede el valor crítico. En general, encontramos que la TRL funciona mejor que la prueba G basada en cualquiera de los métodos de llamada de genotipo (Figura 5). Por ejemplo, a una tasa de falsos positivos del 5% y con una profundidad de secuenciación de 5×, la potencia para detectar un locus de enfermedad con un MAF del 1% y un riesgo relativo (RR) de 2 es del 51% con el LRT, pero la potencia disminuye al 33% utilizando el método de llamada sin filtrado y al 34% utilizando el método de llamada con filtrado. En particular, a baja profundidad, la prueba G aplicada a los llamados genotipos con filtrado funciona muy mal (columna de la izquierda en la Figura 5). Si comparamos la potencia de la LRT con la prueba de tendencia de Armitage utilizando los llamados genotipos, encontramos que la LRT también tiene una potencia más alta que la prueba de tendencia de Armitage (archivo adicional 3). Esto sugiere que si uno desea usar genotipos llamados, filtrarlos basados en la confianza de la llamada puede resultar en una pérdida de potencia.

Figura 5
figura 5

Receiver operating characteristic (ROC) de las cuatro pruebas de asociación. Para la definición de las cuatro estadísticas, véase el pie de foto de la Figura 4. Suponiendo 500 casos y 500 controles, se simuló un conjunto de 20.000 sitios bajo el null y bajo la alternativa a profundidades de secuenciación individuales de 2×, 5× y 10× (tres columnas). A cada tasa de falsos positivos (eje x), el valor crítico correspondiente se calculó utilizando la distribución nula empírica. La tasa positiva verdadera (potencia; eje y) se obtuvo computando la fracción de sitios causales con estadísticas de prueba que exceden el valor crítico.

Aplicación a datos reales

Analizamos 200 exomas de controles para un estudio de asociación de enfermedades que se han secuenciado utilizando la tecnología Illumina a una profundidad por individuo de 8× . Utilizamos las probabilidades de genotipo generadas por el programa» SOAPsnp » para nuestra inferencia. Para obtener más información, consulte Métodos.

En primer lugar, exploramos la precisión de las estimaciones de los MAF a partir de datos de secuenciación de próxima generación para 50 SNP comparándolos con los MAF estimados a partir de datos de genotipo de secuenciación. Tanto las estimaciones que utilizan el método ML como el método de llamada de genotipos sin filtrado están altamente correlacionadas con las estimaciones hechas a partir de los datos de genotipos de secuencia (es decir, una pequeña diferencia estandarizada entre las dos estimaciones en la Figura 6). Sin embargo, las estimaciones basadas en llamadas de genotipos con filtrado muestran una correspondencia pobre con las frecuencias estimadas a partir de los datos de genotipos de secuenciación, especialmente cuando la profundidad de secuenciación es baja. Curiosamente, hay un SNP donde el MAF estimado a partir de los datos de resecuenciación es muy diferente de la estimación obtenida a partir de los datos del genotipo de secuenciación, a pesar de que la profundidad de secuenciación es muy alta (14×). Específicamente, la MAF estimada a partir de los datos del genotipo del secuénomo es de 22,5%, pero es de 17,2% cuando se estima utilizando el abordaje de ML. El examen individual muestra que en muchos individuos, el genotipo altamente respaldado basado en los datos de secuenciación difiere de los genotipos de secuenciación. Dado que este SNP está cubierto por muchas lecturas en estos individuos y que las bases de lectura observadas tienen puntajes de alta calidad (>Q20), es probable que la diferencia se deba a errores de genotipado de secuencia. Tenga en cuenta que hay un par de SNP en los que los MAF estimados del enfoque de llamada de genotipo sin filtrado parecen corresponder mejor a los MAF estimados del genotipado de secuencia que las estimaciones del enfoque de ML. Por ejemplo, en un SNP, el MAF estimado es del 25,7% a partir de los datos del genotipo del Secuénomo25.9% del método de llamada al genotipo sin filtrado y 27,2% del método ML. Sin embargo, la inspección individual revela que hay algunos individuos para los que el genotipo llamado de los datos de secuenciación difiere del genotipo del secuenoma. En estos casos, los errores en los llamados genotipos se cancelaron, dando la apariencia de una mejor correspondencia con los datos del genotipo Secuencial. Por lo tanto, para estos SNPs, es difícil saber qué método funciona mejor.

Figura 6
figura 6

Estimaciones de la frecuencia de los alelos calculadas a partir de 200 individuos utilizando datos de secuenciación de próxima generación vs.datos de genotipo de secuenciación. En cada sitio, solo se utilizaron individuos que tenían datos de genotipo de secuenciación y datos de secuenciación para estimar la frecuencia de los alelos. Para los datos de secuenciación, se obtuvieron estimaciones de MAF utilizando tres métodos diferentes (Llamar NF, Llamar F y ML). La diferencia estandarizada para cada estimación se calculó como , donde y son los MAF estimados a partir de los datos de secuenciación y los datos de genotipo de secuenciación, respectivamente, y n es el número de individuos utilizados para la estimación. Cada sitio se clasifica en uno de cuatro contenedores según la profundidad de cobertura individual media( color): menos de 4×, más de 4× pero menos de 8×, más de 8× pero menos de 16× y más de 16×.

A continuación, examinamos la distribución de los MAF calculados utilizando varios enfoques en un rango de profundidades de secuenciación a partir de nuestros datos de secuenciación de exomas de próxima generación (Figura 7). Descartamos SNPs con un MAF estimado <2%, ya que es difícil distinguir estos SNPs de muy baja frecuencia de los errores de secuenciación en este conjunto de datos. Eliminamos además los sitios en los que había una diferencia significativa (valor de p inferior a 10-5 utilizando una prueba de suma de rangos ) en la puntuación de calidad de las bases de lectura entre los alelos menor y mayor. Es probable que estos sitios sean SNP artificiales que pueden ocurrir debido a un mapeo incorrecto o sesgos desconocidos introducidos durante el procedimiento experimental. Luego clasificamos cada sitio en contenedores basados en la profundidad de cobertura. El número de SNPs en cada bandeja se muestra en la Tabla 1. Cuando la profundidad media es inferior a 9×, las distribuciones de MAF estimados basadas en métodos de llamada a genotipos son muy diferentes de la basada en el método ML. En concreto, los enfoques de llamada de genotipos dan lugar a un gran exceso de SNP de baja frecuencia (MAF entre el 2% y el 4%). Este patrón refleja lo que se vio en nuestros estudios de simulación (Figura 3). Además, para los métodos de llamada de genotipos, la distribución de frecuencia de los alelos cambia drásticamente a medida que cambia la profundidad de secuenciación. Por lo tanto, como se discutió anteriormente, cuando la profundidad no es muy alta, es probable que los métodos de llamada de genotipado incluyan muchos SNPs falsos que son errores de secuenciación. Estos errores aparecen como un exceso de SNPs de baja frecuencia en la distribución de frecuencias. La distribución basada en el método ML es más estable a través de profundidades, pero todavía hay un exceso de SNPs con baja frecuencia alélica con una profundidad inferior a 9×en comparación con la proporción de SNPs de baja frecuencia a mayores profundidades.

Figura 7
figure7

la Distribución de la frecuencia del alelo menor estimarse a partir de los exomas de 200 individuos secuenciados. Para cada sitio, la frecuencia de alelos menores se estimó utilizando cuatro métodos diferentes: (1) el método ML con alelo menor desconocido, (2) el método ML con un alelo menor conocido o fijo, (3) genotipos de llamada sin filtrado (Llamada NF), y (4) genotipos de llamada con filtrado (Llamada F). Cada sitio se clasifica en contenedores según la profundidad de cobertura. Además, en cada histograma, no se consideran los sitios con un MAF estimado inferior al 2%. Para el número de SNP que se utilizaron para este análisis, véase la Tabla 1.

Tabla 1 Número de SNPs con un MAF estimado superior al 2% utilizando un método particular (fila) dentro de cada bin (columna) definido por la profundidad de secuencia promedio entre individuos.

Finalmente, utilizamos estos datos de resecuenciación de exomas para simular un estudio de asociación de casos y controles. Para examinar la distribución de las estadísticas de la prueba de asociación bajo la hipótesis nula, asignamos aleatoriamente 100 individuos a un grupo de casos y los otros 100 al grupo de control. Para todos los SNP en el cromosoma 2 con estimaciones de MAF > 2% (basado en el método desconocido de ML de alelos menores), probamos las diferencias de frecuencia de alelos entre casos y controles computando el estadístico G utilizando genotipos llamados con y sin filtrado, así como el estadístico LRT. La figura 8 muestra los gráficos QQ comparando las distribuciones de las estadísticas de prueba con la distribución estándar χ2 (1). Como se ha visto en estudios de simulación, la distribución nula del estadístico G calculada al llamar genotipos sin filtrar se aparta sustancialmente de la distribución χ2 (1). Sin embargo, la distribución nula de la estadística de TRL sigue de cerca la distribución χ2(1). El factor de inflación es de 1,01, lo que implica que la estadística de TRL funciona bien cuando se aplica a datos reales.

Figura 8
figura 8

QQ-gráficas que comparan las estadísticas de las pruebas de asociación para las diferencias de frecuencia de alelos entre 100 casos y 100 controles con una distribución χ2(1). Los fenotipos se asignaron aleatoriamente a individuos en el conjunto de datos de resecuenciación del exoma, de modo que hay 100 casos y 100 controles. Para cada sitio, se calcularon tres estadísticas: el estadístico G usando genotipos llamados sin filtrado (Llamada NF), el estadístico G-estático usando genotipos llamados con filtrado (Llamada F), y el estadístico LRT. Para minimizar la inclusión de SNP falsos, se descartan los sitios con estimaciones de MAF de ML inferiores al 2%. Para fines de visualización, se muestran los resultados de sitios en el cromosoma 2. Tenga en cuenta que el factor de inflación se muestra en la esquina superior izquierda de cada gráfico QQ.

Related Posts

Deja una respuesta

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