Para todas las preguntas debe mostrar las evidencias procedimentales para cada caso mostrando los códigos usados, la salida de los códigos, los gráficos y las respectivas interpretaciones de los valores.
De acuerdo con los datos de pájaros usados para la primera tarea del laboratorio, conteste las siguientes preguntas y ejercicios:
# Cargar librerías
library(poppr)
library(hierfstat)
library(pegas)
library(tidyr)
library(mmod)
library(data.table)
library(MASS)
library(ggplot2)
library(vcfR)
library(adegenet)
# Cargar base de datos
getwd()
## [1] "/Users/waiwai/Library/Mobile Documents/com~apple~CloudDocs/Universidad/Genetica–poblaciones/Laboratorio-GP-R/02_r_scripts"
pajaros <- read.vcfR("../00_raw_data/datos-Lab01-2/Pajaros.vcf") # Path relativo
## Scanning file to determine attributes.
## File attributes:
## meta lines: 10
## header_line: 11
## variant count: 1925
## column count: 49
##
Meta line 10 read in.
## All meta lines processed.
## gt matrix initialized.
## Character matrix gt created.
## Character matrix gt rows: 1925
## Character matrix gt cols: 49
## skip: 0
## nrows: 1925
## row_num: 0
##
Processed variant 1000
Processed variant: 1925
## All variants processed
pajaros
## ***** Object of Class vcfR *****
## 40 samples
## 594 CHROMs
## 1,925 variants
## Object size: 1.6 Mb
## 0 percent missing data
## ***** ***** *****
# Crear objeto genind
gen.pajaros <- vcfR2genind(pajaros)
# Incorporar estratos
estratos.pajaros <-read.delim("../00_raw_data/datos-Lab01-2/Sporo_pop.txt") #Tabla de la población para hacer estratos
estratos.pajaros
## ID poblacion
## 1 Sc-coDO-727 Sarapiqui
## 2 Sc-coDO-743 Sarapiqui
## 3 Sc-coDO-744 Sarapiqui
## 4 Sc-coUCR042 Sarapiqui
## 5 Sc-coUCR077 Sarapiqui
## 6 Sc-coUM-014 Sarapiqui
## 7 Sc-coUM-015 Sarapiqui
## 8 Sc-coUM-016 Sarapiqui
## 9 Sc-coUM-017 Sarapiqui
## 10 Sc-coUM-018 Sarapiqui
## 11 Sc-coDO-766 Cahuita
## 12 Sc-coDO-767 Cahuita
## 13 Sc-coDO-768 Cahuita
## 14 Sc-coDO-769 Cahuita
## 15 Sc-coUM-046 Cahuita
## 16 Sc-coUM-048 Cahuita
## 17 Sc-coUM-049 Cahuita
## 18 Sc-coUM-050 Cahuita
## 19 Sc-coUM-051 Cahuita
## 20 Sc-coUM-053 Cahuita
## 21 Sc-hoDO-756 Golfito
## 22 Sc-hoDO-757 Golfito
## 23 Sc-hoDO-758 Golfito
## 24 Sc-hoDO-759 Golfito
## 25 Sc-hoLS-115 Golfito
## 26 Sc-hoLS-120 Golfito
## 27 Sc-hoLS-121 Golfito
## 28 Sc-hoUCR586 Golfito
## 29 Sc-hoUCR587 Golfito
## 30 Sc-hoUM-037 Golfito
## 31 Sc-hoDO-760 Quepos
## 32 Sc-hoDO-761 Quepos
## 33 Sc-hoDO-762 Quepos
## 34 Sc-hoDO-763 Quepos
## 35 Sc-hoDO-764 Quepos
## 36 Sc-hoUCR102 Quepos
## 37 Sc-hoUM-038 Quepos
## 38 Sc-hoUM-039 Quepos
## 39 Sc-hoUM-040 Quepos
## 40 Sc-hoUM-043 Quepos
strata(gen.pajaros) <- estratos.pajaros
strata(gen.pajaros)
## ID poblacion
## 1 Sc-coDO-727 Sarapiqui
## 2 Sc-coDO-743 Sarapiqui
## 3 Sc-coDO-744 Sarapiqui
## 4 Sc-coUCR042 Sarapiqui
## 5 Sc-coUCR077 Sarapiqui
## 6 Sc-coUM-014 Sarapiqui
## 7 Sc-coUM-015 Sarapiqui
## 8 Sc-coUM-016 Sarapiqui
## 9 Sc-coUM-017 Sarapiqui
## 10 Sc-coUM-018 Sarapiqui
## 11 Sc-coDO-766 Cahuita
## 12 Sc-coDO-767 Cahuita
## 13 Sc-coDO-768 Cahuita
## 14 Sc-coDO-769 Cahuita
## 15 Sc-coUM-046 Cahuita
## 16 Sc-coUM-048 Cahuita
## 17 Sc-coUM-049 Cahuita
## 18 Sc-coUM-050 Cahuita
## 19 Sc-coUM-051 Cahuita
## 20 Sc-coUM-053 Cahuita
## 21 Sc-hoDO-756 Golfito
## 22 Sc-hoDO-757 Golfito
## 23 Sc-hoDO-758 Golfito
## 24 Sc-hoDO-759 Golfito
## 25 Sc-hoLS-115 Golfito
## 26 Sc-hoLS-120 Golfito
## 27 Sc-hoLS-121 Golfito
## 28 Sc-hoUCR586 Golfito
## 29 Sc-hoUCR587 Golfito
## 30 Sc-hoUM-037 Golfito
## 31 Sc-hoDO-760 Quepos
## 32 Sc-hoDO-761 Quepos
## 33 Sc-hoDO-762 Quepos
## 34 Sc-hoDO-763 Quepos
## 35 Sc-hoDO-764 Quepos
## 36 Sc-hoUCR102 Quepos
## 37 Sc-hoUM-038 Quepos
## 38 Sc-hoUM-039 Quepos
## 39 Sc-hoUM-040 Quepos
## 40 Sc-hoUM-043 Quepos
setPop(gen.pajaros) <- ~poblacion
# Valores perdidos
perdidos <- info_table(gen.pajaros, type = "missing")
## No Missing Data Found!
perdidos #No hay datos perdidos
## NULL
Las estadísticas F indican una estructura poblacional baja. Estas estimativas permiten determinar diferencias en las frecuencias alélicas entre poblaciones en marcadores bialelicos.
El coeficiente de endogamia nos indica de esta población muestra un exceso de heterocigotos promedio de individuos entre las subpoblaciones (Fis < 0; Fis = -0.3347), esto se traduce en apareamientos aleatorios en las subpoblaciones con una probabilidad no significativa de obtener alelos identicos por descendencia (IBD) en la metapoblación. Mientas que el índice de fijación nos indica el déficit de heterocigotos por deriva. Esta población tiene un efecto Wahlund casi nulo, menor al 1% (Fst < 0.2; Fst = 0.008), porque tiene un exceso de heterocigotos que nos indican baja estructura poblacional. Los valores negativos para las estimativas Gst indican un error, ya que los valores estimados se encuentran entre 0-1. Esto podría deberse a la escogencia de marcadores con alta tasa mutacional, o un error de genotipeo.Gst se implementa para marcadores multialélicos, en este caso los loci por individuo tiene dos alelos por ser SNP’s. No se recomienda su uso.
Se obtienen un Gst globales negativo (Gst_Nei = -0.0132; Gst_Hedrick = -0.0258), y los datos no permiten el cálculo de significancia por medio de bootstraps. Por tanto, se puede interpretar como la ausencia de diferencias en las frecuencias alélicas entre las subpoblaciones, sin poder asegurar por medio de límites de confianza que Gst de Nei, ni Gst de Headrick sean cero. Asumiendo que los Gst sean cero, el resultado sugiere que existe migración y flujo entre las subpoblaciones que reduce los efectos de la deriva sobre la diversidad de las poblaciones.
En una publicación, indicaría las estadísticas F, ya que estos son màs apropiados para analizar marcadores bialélicos y han sido ampliamente usados en la literatura para determinar y comparar distancias genéticas, presentado entonces una ventaja por familiaridad. Mientras que las estadísticas Gst son susceptibles a tasas de mutación altas, marcadores multialélicos y son más complicadas de interpretar. Gst de Hedrick es una medida estandarizada que soluciona algunes errores que puede producir la Gst de Nei, sin embargo, aún no es confiable para esta clase de marcador.
F.pajaros <- Fst(as.loci(gen.pajaros)) #"as.loci" de pegas transforma archivos genind
head(F.pajaros, 6)
## Fit Fst Fis
## SCONTIG0_55544 -0.22615804 -0.02633969 -0.19469027
## SCONTIG0_652706 -0.10655738 0.03005464 -0.14084507
## SCONTIG0_1163450 -0.10599078 -0.02918587 -0.07462687
## SCONTIG0_1475736 -0.04347826 -0.01449275 -0.02857143
## SCONTIG1_857604 -0.47727273 -0.01641414 -0.45341615
## SCONTIG1_912009 -0.09589041 0.01065449 -0.10769231
colMeans(F.pajaros) #promedios por columna
## Fit Fst Fis
## -0.324314394 0.008053526 -0.334697792
# Fst pareado por población
pf <- pairwise.fst.dosage(gen.pajaros, pop = gen.pajaros$strata$poblacion) #lo hace por población
pf
## Sarapiqui Cahuita Golfito Quepos
## Sarapiqui NA 0.016127369 0.004510127 0.002776518
## Cahuita 0.016127369 NA 0.006109447 0.005472122
## Golfito 0.004510127 0.006109447 NA 0.004033498
## Quepos 0.002776518 0.005472122 0.004033498 NA
Gst_Nei.pajaros <- Gst_Nei(gen.pajaros)
Gst_Nei.pajaros$global
## [1] -0.01317034
Gst_Hedrick.pajaros <- Gst_Hedrick(gen.pajaros)
head(Gst_Hedrick.pajaros$per.locus,6)
## SCONTIG0_55544 SCONTIG0_652706 SCONTIG0_1163450 SCONTIG0_1475736
## -0.053970904 0.028533742 -0.040863198 -0.017746229
## SCONTIG1_857604 SCONTIG1_912009
## -0.075022065 0.006166586
Gst_Hedrick.pajaros$global
## [1] -0.02580427
Noo hay estructura génetica significativa (Fis = -0.3347; Fst = 0.008;Gst_Nei = -0.0132; Gst_Hedrick = -0.0258), interpretado a partir de valores de endogamia, y estructura génica negativos. Sin embargo, desconocemos si en realidad incluyen o no 0. Puesto a que no podemos determinar su significancia a partir de bootstraps, esto simplemente se asume a partir de la interpretación de los datos.
Gst_p.pajaros <- pairwise_Gst_Nei(gen.pajaros)
Gst_p.pajaros
## Sarapiqui Cahuita Golfito
## Cahuita -0.004007168
## Golfito -0.009829955 -0.008989542
## Quepos -0.010582478 -0.009192254 -0.009869710
Determine si existe aislamiento por distancia (IBD) entre las diferentes poblaciones de muestreo mediante el uso de la prueba de Mantel.
¿Qué se puede inferir del gráfico? Respalde su respuesta interpretando el valor de correlación y el valor p estadístico de la prueba de Mantel.
dgeo <- read.csv("../00_raw_data/datos-Lab01-2/dist_pajaros.csv", row.names=1, header = TRUE, sep = ",")
dgeo
## x y
## Sarapiqui -84.00746 10.430218
## Cahuita -82.82082 9.714834
## Golfito -83.11154 8.602726
## Quepos -84.14660 9.432652
#Calcular distancias euclideanas geograficas
dist.pajaros <- dist(dgeo)
dist.pajaros
## Sarapiqui Cahuita Golfito
## Cahuita 1.385601
## Golfito 2.035289 1.149479
## Quepos 1.007222 1.355476 1.326695
#Calcular distancias genéticas
dist.gen <- pairwise_Gst_Nei(gen.pajaros)
dist.gen
## Sarapiqui Cahuita Golfito
## Cahuita -0.004007168
## Golfito -0.009829955 -0.008989542
## Quepos -0.010582478 -0.009192254 -0.009869710
# Prueba de Mantel permite comparar entre matrices
ibd <- mantel.randtest(dist.pajaros, dist.gen, nrepet=999)
ibd
## Monte-Carlo test
## Call: mantel.randtest(m1 = dist.pajaros, m2 = dist.gen, nrepet = 999)
##
## Observation: 0.03027893
##
## Based on 999 replicates
## Simulated p-value: 0.406
## Alternative hypothesis: greater
##
## Std.Obs Expectation Variance
## -0.007579195 0.033768255 0.211951968
La prueba de Monte-Carlo permite realizar una comparación entre las distancias euclideanas geográficas entre poblaciones y sus distancias genéticas. Se obtiene una correlación baja entre distancia geográfica y genética ( r_m < 0.3; r_m = 0.03 p = 0.36 ). Esto quiere decir que un 3% de la variación de los datos son explicados por la relación entre distancias geográficas y genéticas, sin embargo, no son diferente a cero por lo que no existe correlación. No hay evidencia de impedimentos en el flujo génico entre las poblaciones de pájaros.
df<-data.frame(dist_geo = as.vector(dist.pajaros), # Para linearizar las distancias en una base de datos
dist_gen = as.vector(dist.gen))
ggplot(df, aes(dist_geo, dist_gen))+ #definimos x y y
geom_point()+ #ponemos puntos encima
geom_smooth(method = "lm", se = F, color = "red", linetype = "dashed")+ #ponemos una regresion lineal, sin error estandar
theme_classic() + xlab("Distancia geográfica") + ylab("Gst pareado") +
theme(axis.text.y = element_text(angle = 90, vjust = 0.5, hjust = 1))
## `geom_smooth()` using formula = 'y ~ x'
Fig. 1 Relación lineal entre Gst de Nei pareado de las poblaciones de
pájaros, y la distancia geográfica con transformación euclideana.
Muestra una pendiente muy baja, lo cual indica que no existe estructura
genética entre las poblaciones.
A partir del cálculo pareado de Gst, se observa una mayor estructura génica entre Quepos y Sarapiquí (Gst_Nei_Q-S = -0.010), que indica menor flujo génico; y menor Gst entre Cahuita y Sarapiquí (Gst_Nei_C-S =-0.004) que indica que hay más migrantes entre estos sitios. Los demás sitios presentan una estructura similar, entre 0.008-0.009.
Sin embargo, puesto a que todos son negativos y no se pueden calcular límites de confianza, se asume que no existe del todo estructura espacial entre las poblaciones pareadas.
arbol.pob <- pairwise_Gst_Nei(gen.pajaros)%>%
nj() %>% # calcula el árbol nj
ladderize() # organiza las ramas por población
plot(arbol.pob, type ="radial", font= 1, cex = 0.5) #realiza el gráfico
Fig 2. Dendrograma de relación entre distancias geográficas y distancias
genéticas en términos de Gst de Nei.
Del dendrograma se puede concluir que existe mayor flujo y similitud entre las poblaciones de Quepos y Sarapiquí, con respecto a Golfito y Cahuita como grupo. Existe mayor similitud también entre Golfito y Cahuita con respecto a las demás poblaciones, sin embargo la diferencia entre Golfito y Cahuita es mayor que la que hay entre el otro grupo de poblaciones. Si se compara este dendrograma con el Gst de Nei pareado, esto se cumple.