Post on 09-May-2018
Frecuencia de caras al lanzar una moneda sucesivamente
Ley de los grandes número
Demostrar mediante simulación la ley de los grandes número Nota: la demostración puede hacerse simulando el lanzamiento de un dado o moneda, el resultado de un cruce mendeliano,…
> prob = c(0.5,0.5) #lanzamiento una moneda misma probabilidad > moneda = c(1,0) #1 cara y 0 cruz (valores numéricos para obtener probabilidad directamente) > nsample = c(1:200) #los datos muestreados se incrementa en cada nuevo muestreo > nreplica = 200 #número muestras simuladas > vector_media = c(1:nreplica) > for (i in 1: nreplica){vector_media[i]=mean(sample(moneda,nsample[i],prob,replace=T))} #bucle vector probabilidad > vector_media > plot(vector_media, main="Ley grandes números", xlab="Número tiradas", ylab="Proporción Caras")
Solución: Ley grandes números demostrada mediante simulación
de la frecuencia de caras al lanzar una moneda sucesivamente
Bucle (Loop) for (i in start:finish) Ejecuta función >y =c(1:10) >for (i in 1:10){ y[i]= sqrt(i)}
THE DATA PIPELINE Probability statistical models
Probability distributions
> help(distribution) # también ?distribution
Probability distributions in R
Los modelos son caricaturas de la realidad
MODELOS PROBABILÍSTICOS
http://www.johndcook.com/blog/distribution_chart/
Averigua la distribución Indique qué clase de distribución siguen probablemente las siguientes variables. Elija siempre la solución más precisa.
1: ¿Producción lechera, en litros/vaca y día?
2: ¿Número de ejemplares de cada especie en una muestra de 20 insectos pertenecientes a una comunidad de 4 especies distintas?
3: ¿Frecuencia de machos en una muestra de 50 individuos de una especie bisexual?
4: ¿(Altura de un individuo - promedio en la población)/desviación típica?
7: ¿Número de hijos por familia en una población humana?
5: El número de entrecruzamientos por cromosoma y meiosis
Distribución normal
2)(2
1
2
1)(
x
exf
# Poisson distribution > ppois(16, lambda=12) # cumulative lower tail > ppois(16, lambda=12, lower=FALSE) # upper tail > pois = rpois(100, lambda=3) # muestra aleatoria tamaño 100 > qpois(.95,3) #value n 95% cumulative lower tail
# Continuous Uniform Distribution
> runif(10, min=1, max=3) # 10 valores entre 1 y 3
# Exponential Distribution
> pexp(2, rate=1/3) #acumulative lower tail > rexp(12, rate=1/3) #muestra aleatoria tamaño 12
Distribuciones discretas
# binomial > dbinom(4, size=12, prob=0.2) > binomial <- rbinom(100, 20, 0.2) # n resultados favorables p = 0.2 en muestras de m=20 y 100 simulaciones
La distribución normal
2)(2
1
2
1)(
x
exf
Dos parámetros determinan la forma de la distribución:
•la media, y
•la desviación típica (raíz cuadrada de
la varianza 2 = [(x-)2]/n)
•A y B difieren en sus medias (4 y 8) •B y C difieren en sus desviaciones típicas (1 y 0,5)
Función de densidad
Normal Distribution
# density normal distribution > dnorm (-1.1) #valor f(x) de N(0,1) #plot density function of normal distribution N(0,1) > x <- seq(-4, 4, 0.1) > plot(x,dnorm(x), type = "l") #plot density function of normal distribution N(2,.25) > x <- seq(0, 4, 0.1) > plot(x,dnorm(x,2,0.25), type = "l") # N(2,0.25) #Cumulative probability > pnorm(84, mean=72, sd=15.2, lower.tail=FALSE) #cumulative N(72,15.2) #Quantile value > qnorm(.975) #z value 97.5% cumulative lower tail # Sampling normal distribution > muestranorm = rnorm(100) #sample n=100 of z values > mean(muestranorm) > sd(muestranorm) >normalmea2s3 = rnorm(300,2,3) #sample n=300 of N(2,3) >mean(normalmea2s3)
The R Project for Statistical Computing
2)(2
1
2
1)(
x
exf
Q-Q Plot
> x<-c(4.75, 3.4, 1.8, 2.9, 2.2, 2.4, 5.8, 2.6, 2.4, 5.25) > qqnorm(x) > qqline(x)
Q-Q Plot
Line of quantiles of the data against normal quantiles
Tema 16: Genética cuantitativa 73
Muestreo de una
distribución normal
Observed data point
P(x = 2,5) ? P(x ≥ 2,5); μ = 3; σ = 2 P(x ≥ 185); μ = 170; σ = 12 P(-1 ≤ z ≤ 2,5); μ = 0; σ = 1
http://www.mathsisfun.com/data/standard-normal-distribution-table.html
Standard score =>
2)(2
1
2
1)(
x
exf2
2
1
2
1)(
x
ezf
THE DATA PIPELINE Summary statistics
•Estadístico (cantidad calculada de la muestra: media, varianza muestral) ( , s2)
•Parámetro ( , σ2 ) µ
•Distribución muestral de un estadístico
•Distribución de la variable
•Distribución del estadístico muestral
Inferencia estadística
Propiedades de los estimadores
•Sesgo
•Eficiencia
•Consistencia
THE DATA PIPELINE Summary statistics
•Estadístico (cantidad calculada de la muestra: media, varianza muestral) ( , s2)
•Parámetro ( , σ2 ) µ
•Distribución muestral de un estadístico
•Distribución de la variable
•Distribución del estadístico muestral
s estimador de σ
Distribución muestral de un estadístico
Distribución muestral de un estadístico
Promedio Promedio Promedio
The distribution of sample means from most distributions will be approximately normally distributed
The distribution of sample means from most distributions will be approximately normally distributed
http://www.nature.com/nmeth/journal/v10/n9/fig_tab/nmeth.2613_F3.html
# Sampling distribution of the statistics mean for a normal distributed variable > n_muestra = 100 > n_replicas = 1000 > vector_medias_muestras = seq(1,n_replicas) #vector_medias_muestra por réplica #bucle que calcula las n_replicas medias de cada muestra simulada de tamaño n_muestra # de una distribución normal > for (i in 1: n_replicas){vector_medias_muestras[i]=mean(rnorm(n_muestra))} > vector_medias_muestras > hist(vector_medias_muestras, nclass = 12) > mean(vector_medias_muestras) > sd(vector_medias_muestras) > 1/sd(vector_medias_muestras) # Sampling distribution of the mean statistics for an uniform distribution > for (i in 1: n_replicas){vector_medias_muestras[i]=mean(runif(n_muestra , min=1, max=3))} > hist(vector_medias_muestras) # Sampling distribution of the mean statistics for an exponential distribution > for (i in 1: n_replicas){vector_medias_muestras[i]=mean(rexp(100))} > hist(vector_medias_muestras)
Simulación de la distribución muestral del estadístico media en tres distribuciones de probabilidad: normal, uniforme y exponencial
Error del muestreo
Distribución muestral de la media muestral
D Distribución x Distribución
Error estándar de la media (Estándar error of mean)
85
Teorema Central del Límite
86
Teorema Central del Límite
La distribución aproximadamente normal de muchas variables muestrales puede fundamentarse en el teorema central del límite, que dice que la distribución de la suma aleatoria de muchos pequeños efectos se aproxima asintóticamente a una distribución normal, “a singularly beautiful law” (F. Galton).
87
Galton’s apparatus (Galton’s Quincunx)
Demostración intuitiva del Teorema Central
del Límite
Ver animación aquí
La posición que ocupa cada bola en el eje X resulta de la suma de pequeños efectos aleatorios +/- y su distribución se aproxima asintóticamente a una distribución normal.
88
Ver animación aquí
Galton’s apparatus (Galton’s Quincunx) Demostración intuitiva
del Teorema Central del Límite
La posición que ocupa cada bola en el eje X resulta de la suma de pequeños efectos aleatorios +/- y su distribución se aproxima asintóticamente a una distribución normal.
89
Significado biológico: cualquier variable influida por la acumulación de multitud de causas con un efecto pequeño y aleatorio tenderán a distribuirse normalmente
Significado estadístico: Muchos de los estimadores que se usan para hacer inferencias de parámetros de la población son sumas o promedios de medidas de la muestra, por lo que podemos usar la distribución normal para describir su conducta y efectuar inferencias
Teorema Central del Límite
90
Leyes universales y TLC (teorema límite central) Una ley macroscópica universal es aquella que gobierna el comportamiento colectivo de un sistema con independencia de su estructura microscópica. La ley de los grandes números o el teorema del límite central son casos de leyes universales cuyo significado se entiende. La causa matemática que subyace a la universalidad de otras leyes se desconoce
Teorema Central del Límite
THE DATA PIPELINE Inferencia
Estimación puntual e intervalo de confianza de un parámetro poblacional
http://2012books.lardbucket.org/books/beginning-statistics/s11-estimation.html#fwk-shafer-ch07_s01_f02
http://2012books.lardbucket.org/books/beginning-statistics/s11-estimation.html#fwk-shafer-ch07_s01_f02
#Zalpha/2 >alpha = 0.05 >qnorm(1-alpha/2)
Compara intervalos de confianza de distribución z para n=100 alpha=0.1 alpha=0.05 alpha=0.01 alpha=0.001
Estimación por intervalos de
confianza
Computer Simulation of 40 confidence Intervals of a mean where α = 0.05
n-1 dfs
n ≥ 30 t -> z
n < 30
x <- seq(-8,8,by=.1)
plot(x,dnorm(x),type='l',ylab="",main="df=2") lines(x,dt(x,df=2),lty=2)
Distribución t de Student
n-1 dfs
n ≥ 30 t -> z
n < 30
#talpha/2 >alpha = 0.05 >n = 5 >qt(1-alpha/2 , df=n−1)
Distribución t de Student
# Intervalo de confianza #Intervalo de confianza del 95% para la variable GC% #Se debe calcular n, media (xbar), desv típica (s), error estándar media (sem) n = sapply(GC,length) xbar = sapply(GC,mean) s = sapply(GC,sd) # sample standard deviation sem = s/sqrt(n); # standard error of mean E = qnorm(.975)∗sem; # margin of error - If t student then E = qt(.975, df=n−1)∗sem; IC = xbar + c(-E,E)
Datos Tabla resumen genoma humano
Loc Type Name RefSeq Size (Mb) GC% Protein rRNA tRNA Other RNA Gene Pseudogene
Nuc Chr 1 NC_000001.
11 248.96 42.3 9,684 17 89 3,889 4,778 1,050
Nuc Chr 2 NC_000002.
12 242.19 40.3 6,720 - 7 3,286 3,600 853
Nuc Chr 3 NC_000003.
12 198.3 39.7 5,616 - 4 2,403 2,780 668
Nuc Chr 4 NC_000004.
12 190.22 38.3 3,817 - 1 2,072 2,251 581
Nuc Chr 5 NC_000005.
10 181.54 39.5 3,949 - 17 2,017 2,393 579
Nuc Chr 6 NC_000006.
12 170.81 39.6 4,674 - 143 2,187 2,828 685
Un - . 133.2 44.6 8,128 4 184 3,042 5,260 1,424
•Estima intervalo confianza media alpha=0.05 siguiente conjunto de datos data = c(0.467, 0.645, 0.868, 0.472, 0.844, 0.879, 0.405, 0.604, 0.787, 0.449, 0.772, 0.780)
•Comparar valores con para alpha =0.05 y n=2,4,6,10,20,30
•El siguiente enlace se encuentran datos de la redundancia de los genomas secuenciados por el proyecto DGRP en D. melanogaster -> Lines_coverage_DGRP Calcula el intervalo de confianza del 95 para la variable Covered_Bases Debes calcular s # sample standard deviation sem = s/sqrt(n); # standard error of mean E = qnorm(.975)∗sem; # margin of error - If t student then E = qt(.975, df=n−1)∗sem; xbar + c(-E,E)
THE DATA PIPELINE Inferencia
Contraste de hipótesis / prueba de hipótesis
Elements of a statistical test
•Null hypothesis: H0
•Alternative hypothesis: H1
•Region of rejection / Critical region
•p-value: the value p of the statistic used to test the null hypothesis
•α: level of significance
Observed data point
Observed data point
Critical value
α
• Elements of a statistical test
•Null hypothesis: H0
•Alternative hypothesis: H1
•Region of rejection / Critical region
•p-value: the value p of the statistic used to test the null hypothesis
•α: level of significance
Statistical hypothesis testing
Two-tailed hypothesis testing One-tailed hypothesis testing
Statistical hypothesis testing
A B
Observed data point
•Errors in making a decision
•Error type I (α): H0 is rejected even though it is true (false positive)
•Error type II (β): H0 is not rejected even though it is false (false negative)
•Power of a test (1 – β): probability to find an effect when it exists
Nombre Fórmula Notas
Test-z para una muestra
(Población distribuida normal o n > 30) y σ conocida. (z es la distancia desde la media en relación con la desviación estándar de la media). Para distribuciones no normales es posible calcular una proporción mínima de una población que cae dentro de k desviaciones estandar para cualquier k.
Test-z para dos muestras Población normal y observaciones independientes con σ1 y σ2 conocidas
t-test para una muestra
(Población normal o n > 30) y desconocida
t-test para datos apareados
(Población normal de diferencias o n > 30) y desconocida o pequeña muestra de tamaño n < 30
Pruebas de hipótesis (ejemplos)
Más tests en -> http://es.wikipedia.org/wiki/Contraste_de_hipótesis
The null hypothesis is that μ = 15.4. We begin with computing the test statistic. > xbar = 14.6 # sample mean > mu0 = 15.4 # hypothesized value > sigma = 2.5 # population standard deviation > n = 35 # sample size > z = (xbar−mu0)/(sigma/sqrt(n)) > z # test statistic [1] −1.8931 We then compute the critical values at .05 significance level. > alpha = .05 > z.half.alpha = qnorm(1−alpha/2) > c(−z.half.alpha, z.half.alpha) [1] −1.9600 1.9600
Two-Tailed Test of Population Mean with Known Variance
n = 35 sigma = 2.5
El siguiente enlace se encuentran datos de la redundancia de los genomas secuenciados por el proyecto DGRP en D. melanogaster -> Lines_coverage_DGRP Prueba si media
Finches_dataset
Measurements in a Sample of Medium Ground Finches from Daphne Major (data were provided by scientists Peter and Rosemary Grant)
Compara las medias de los distintas variables medidas en Geospiza entre los años 1977 y el resto. ¿Hay algún fenotipo sometido a selección?
Multinomial Goodness of Fit
Cruces dihíbrido de Mendel: F2
315 plantas lisas-amarillas 108 lisas-verdes 101 rugosas-amarillas 32 rugosas-verdes a la F2.
Analizar mediantes una prueba de bondad de ajuste de chi-cuadrado si:
(a) Se ajustan a una razón 9:3:3:1 (b) las proporciones lisas:rugosas son 3:1 (c) les proporciones amarillas:verdes son 3:1
> f2observada = c(315,108,101,32)
> f2probteorica = c(9/16,3/16,3/16,1/16)
> chisq.test(f2observada, p=f2probteorica)
Chi-squared Test of Independence
> B = matrix( c(38, 76, 15, 105, 122 , 13), nrow=2, ncol=3, byrow=TRUE) > chisq.test(B)
CT CC TT
Control 38
(29,5%)
76
(58,9%)
15
(11,6%)
Case 105
(43,8%)
122
(50,8%)
13
(5,4%)
> Mapped_Coverage.lm = lm(Mapped_Coverage ~ Covered_Bases, data=coverage) > coeffs = coefficients(Mapped_Coverage.lm); coeffs > summary (Mapped_Coverage.lm) #F-statistic significance test
Single linear regression
# Regresión lineal #Regresión lineal variables Mapped_Coverage y Covered_Bases y de los datos Lines_coverage_DGRP #Código R
Finches_dataset
anova(lm(Beak_Width_(mm) ~ Last_Year,data= finches_dataset))
# ANOVA
Análsis de la varianza
Decision tree statistical test
Type of Data
Goal Measurement (from Gaussian Population)
Rank, Score, or Measurement (from
Non- Gaussian Population)
Binomial (Two Possible
Outcomes)
Survival Time
Describe one group Mean, SD Median, interquartile range
Proportion Kaplan Meier survival curve
Compare one group to a hypothetical value
One-sample t test Wilcoxon test Chi-square or Binomial test **
Compare two unpaired groups
Unpaired t test Mann-Whitney test Fisher's test (chi-square for large samples)
Log-rank test or Mantel-Haenszel*
Compare two paired groups
Paired t test Wilcoxon test McNemar's test Conditional proportional hazards regression*
Compare three or more unmatched groups
One-way ANOVA Kruskal-Wallis test Chi-square test Cox proportional hazard regression**
Compare three or more matched groups
Repeated-measures ANOVA
Friedman test Cochrane Q** Conditional proportional hazards regression**
Quantify association between two variables
Pearson correlation Spearman correlation Contingency coefficients**
Predict value from another measured variable
Simple linear regression or Nonlinear regression
Nonparametric regression**
Simple logistic regression* Cox proportional hazard regression*
Predict value from several measured or binomial variables
Multiple linear regression* or Multiple nonlinear regression**
Multiple logistic regression*
Cox proportional hazard regression*
http://www.graphpad.com/support/faqid/1790/
Choosing a Statistical Test
How to determine the appropriate statistical test • I find that a systematic, step-by-step approach is the best way to decide how to
analyze biological data. I recommend that you follow these steps: • Specify the biological question you are asking. • Put the question in the form of a biological null hypothesis and alternate hypothesis. • Put the question in the form of a statistical null hypothesis and alternate hypothesis. • Determine which variables are relevant to the question. • Determine what kind of variable each one is. • Design an experiment that controls or randomizes the confounding variables. • Based on the number of variables, the kinds of variables, the expected fit to the
parametric assumptions, and the hypothesis to be tested, choose the best statistical test to use.
• If possible, do a power analysis to determine a good sample size for the experiment. • Do the experiment. • Examine the data to see if it meets the assumptions of the statistical test you chose
(primarily normality and homoscedasticity for tests of measurement variables). If it doesn't, choose a more appropriate test.
• Apply the statistical test you chose, and interpret the results. • Communicate your results effectively, usually with a graph or table
Online Handbook of Biostatistics
http://www.r-tutor.com/elementary-statistics
THE DATA PIPELINE Multiple testing
• Bonferroni correction
• False discovery rate
Multiple testing
http://xkcd.com/882/
120
p-hacking (selective reporting or Inflation bias)
Selectively report analyses that produce significant results
• conducting analyses midway through
experiments to decide whether to continue collecting data
• recording many response variables and deciding which to report postanalysis
• deciding whether to include or drop outliers postanalyses, excluding, combining, or splitting treatment groups postanalysis
• including or excluding covariates postanalysis stopping data exploration if an analysis yields a significant p-value
Head ML, Holman L, Lanfear R, Kahn AT, Jennions MD (2015) The Extent and Consequences of P-Hacking in Science. PLoS Biol 13(3): e1002106. doi:10.1371/journal.pbio.1002106
The effect of p-hacking on the distribution of p-values in the range of significance
Links Analysis tools
•Elementary Statistics with R
•Online Handbook of Biostatistics
•R tutorial. An R Introduction to Statistics
•Statistical Analysis in R
•Statistical calculations
•Statistical errors. Features in Nature 13 february 2014
122 Ritchie MD, Holzinger ER, Li R, Pendergrass SA, Kim D. Methods of integrating data to uncover genotype-phenotype interactions. Nature reviews. Genetics. 2015 Feb;16(2):85-97.
123
•Multi-staged approach: divide data analysis into multiple steps, and signals are enriched with each step of the analysis. Reflect hypothesis variation is hierarchical. Example triangle approach of predicting causal genes in GWAS loci
1. SNPs are associated with the phenotype and filtered based on a genome-wide significance threshold 2. SNPs deemed significant from step 1 are then tested for association with another level of omic data (eQTLs, mQTLs, metabolite QTLs, pQTLs) 3. Omic data used in step 2 are then tested for correlation with the phenotype of interest.
• Meta-dimensional approach: Meta-dimensional analysis combines multiple data types in a simultaneous analysis
Data integration methods