setwd("~/Dropbox/Clases_UPJ/Biología de poblaciones/A_Ecopo/Ciclo 1910/Examen_2")
## Punto 1
# Tasa intrínseca de crecimiento
datos1<-read.table("inva.csv", sep = ",", header=T)
datos1
## año n
## 1 1962 2
## 2 1975 200
## 3 1986 264
## 4 1989 328
## 5 1994 2200
## 6 1999 3955
attach(datos1)
plot (año, n)
datos1.tra<-transform(datos1, n=log(n), año=año-1962)
detach(datos1)
attach(datos1.tra)
plot (año, n)
fit1<-lm(n~año)
abline(fit1)
summary(fit1)
##
## Call:
## lm(formula = n ~ año)
##
## Residuals:
## 1 2 3 4 5 6
## -0.69394 1.46597 -0.32546 -0.67268 0.29003 -0.06393
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.38709 0.75809 1.830 0.14127
## año 0.18810 0.02986 6.299 0.00325 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9052 on 4 degrees of freedom
## Multiple R-squared: 0.9084, Adjusted R-squared: 0.8855
## F-statistic: 39.68 on 1 and 4 DF, p-value: 0.003246
# Tiempo doble
(t2<-log(2)/0.18810)
## [1] 3.684993
# Tiempo necesario para alcanzar 10000 individuos
m.time1<-function(r, N=10000) {(log(N)-1.38709)/r}
rs<-c(0.1881*2,0.1881*3,0.1881*4)
m.time1(rs)
## [1] 20.79546 13.86364 10.39773
### Sección II
## A: Matriz de Lefkovitch
A=matrix(c(0,40,70,60,0.25,0,0,0,0,0.4,0,0,0,0,0.25,0), nr=4, byrow = T)
A
## [,1] [,2] [,3] [,4]
## [1,] 0.00 40.0 70.00 60
## [2,] 0.25 0.0 0.00 0
## [3,] 0.00 0.4 0.00 0
## [4,] 0.00 0.0 0.25 0
N0=matrix(c(400,100,250,75), ncol=1)
N0
## [,1]
## [1,] 400
## [2,] 100
## [3,] 250
## [4,] 75
## B: Proyección a 10 años
years=10
n.projections=matrix(0,nrow=nrow(A), ncol=years+1)
n.projections
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
## [1,] 0 0 0 0 0 0 0 0 0 0 0
## [2,] 0 0 0 0 0 0 0 0 0 0 0
## [3,] 0 0 0 0 0 0 0 0 0 0 0
## [4,] 0 0 0 0 0 0 0 0 0 0 0
n.projections[,1]=N0
n.projections
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
## [1,] 400 0 0 0 0 0 0 0 0 0 0
## [2,] 100 0 0 0 0 0 0 0 0 0 0
## [3,] 250 0 0 0 0 0 0 0 0 0 0
## [4,] 75 0 0 0 0 0 0 0 0 0 0
for (i in 1:years) n.projections[,i+1]=A%*%n.projections[,i]
n.projections
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,] 400 26000.0 10550 263400.0 288100 2746850.00 4740625.0 29880300.0
## [2,] 100 100.0 6500 2637.5 65850 72025.00 686712.5 1185156.2
## [3,] 250 40.0 40 2600.0 1055 26340.00 28810.0 274685.0
## [4,] 75 62.5 10 10.0 650 263.75 6585.0 7202.5
## [,9] [,10] [,11]
## [1,] 67066350.00 336107650.0 886936537.5
## [2,] 7470075.00 16766587.5 84026912.5
## [3,] 474062.50 2988030.0 6706635.0
## [4,] 68671.25 118515.6 747007.5
## C: Gráfica por estados
matplot(0:years, t(n.projections), type = "l", lty = 1:4,col = 1, ylab = "Abundancia", xlab = "Año")
legend("topleft", legend = c("Semilla", "Adulto pequeño", "Adulto mediano", "Adulto grande"),lty = 1:4, col = 1, bty = "n")
## D: lambda que converge en el tiempo, usando ánalysis de valores y vectores "propios" (Eigenanalysis) de la matriz llegamos directamente a la convergencia de lambda
eigs.A=eigen(A)
eigs.A
## eigen() decomposition
## $values
## [1] 3.4832784+0.0000000i -2.7690598+0.0000000i -0.3571093+0.1672943i
## [4] -0.3571093-0.1672943i
##
## $vectors
## [,1] [,2] [,3]
## [1,] 0.9974004611+0i -0.995864517+0i 0.7103908+0.0000000i
## [2,] 0.0715848938+0i 0.089909987+0i -0.4078193-0.1910503i
## [3,] 0.0082204045+0i -0.012987800+0i 0.2923829+0.3509685i
## [4,] 0.0005899905+0i 0.001172582+0i -0.0734620-0.2801157i
## [,4]
## [1,] 0.7103908+0.0000000i
## [2,] -0.4078193+0.1910503i
## [3,] 0.2923829-0.3509685i
## [4,] -0.0734620+0.2801157i
dom.pos=which.max(eigs.A[["values"]])
## Warning in which.max(eigs.A[["values"]]): imaginary parts discarded in
## coercion
L1 <- Re(eigs.A[["values"]][dom.pos])
L1
## [1] 3.483278
# Años necesarios para la convergencia: tasa de crecimiento anual (El lambda general de la población es el valor que converge en el tiempo), lambda=Nt+1/Nt para cada año.
# Primero necesitamos la suma de todos los estados aplicando la función sum a cada columna.
n.totals=(apply(n.projections,2,sum))
n.totals
## [1] 825.0 26202.5 17100.0 268647.5 355655.0
## [6] 2845478.8 5462732.5 31347343.8 75079158.8 355980783.1
## [11] 978417092.5
# Ahora obtenemos cada Lambda dividiendo todos los Nt+1 por cada Nt.
lambdas=n.totals[-1]/n.totals[-(years+1)]
lambdas
## [1] 31.7606061 0.6526095 15.7103801 1.3238724 8.0006713 1.9197938
## [7] 5.7383999 2.3950724 4.7414061 2.7485110
#Los lambda que obtuvimos no convergen en el tiempo entonces necesitamos más unidades de tiempo
t=90
Nt=N0/sum(N0)
L.t=numeric(t)
#Creamos una iteración (for-loop) que reusa Nt para cada paso del tiempo, asegurandonos que tengamos un vector numerico vacio para sostener las salidas
for (i in 1:t) L.t[i] <- {
Nt1 <- A %*% Nt
L <- sum(Nt1)/sum(Nt)
Nt <- Nt1/sum(Nt1)
L
}
L.t # La convergencia se alcanza en el año 84
## [1] 31.7606061 0.6526095 15.7103801 1.3238724 8.0006713 1.9197938
## [7] 5.7383999 2.3950724 4.7414061 2.7485110 4.2235392 2.9979446
## [13] 3.9315583 3.1675477 3.7592892 3.2799713 3.6549170 3.3532406
## [19] 3.5906619 3.4004660 3.5507140 3.4306881 3.5257264 3.4499403
## [25] 3.5100369 3.4621687 3.5001621 3.4699214 3.4939376 3.4748307
## [31] 3.4900103 3.4779371 3.4875310 3.4799019 3.4859652 3.4811442
## [37] 3.4849761 3.4819295 3.4843512 3.4824259 3.4839563 3.4827396
## [43] 3.4837068 3.4829379 3.4835491 3.4830632 3.4834495 3.4831424
## [49] 3.4833865 3.4831925 3.4833467 3.4832241 3.4833216 3.4832441
## [55] 3.4833057 3.4832567 3.4832957 3.4832647 3.4832893 3.4832698
## [61] 3.4832853 3.4832729 3.4832828 3.4832750 3.4832812 3.4832762
## [67] 3.4832802 3.4832770 3.4832795 3.4832775 3.4832791 3.4832779
## [73] 3.4832789 3.4832781 3.4832787 3.4832782 3.4832786 3.4832783
## [79] 3.4832785 3.4832783 3.4832785 3.4832784 3.4832785 3.4832784
## [85] 3.4832784 3.4832784 3.4832784 3.4832784 3.4832784 3.4832784
# Número relativo de individuos en la convergencia, la distribución estable de estados se alcanza cuando la abundancia relativa de los estados es la misma
w <- Re(eigs.A[["vectors"]][, dom.pos])
ssd <- w/sum(w)
round(ssd, 3)
## [1] 0.925 0.066 0.008 0.001
# Lo anterior muestra que si la matriz de projección no cambia en el tiempo, estara compuesta de los porcentajes calculados para cada estado
## E: Las hipótesis deberian plantearse sobre el estado 2 (Adultos Pequeños) y el estado 3 (Adultos Medianos), ya que son los estados con mayor fertilidad.
# Para los Adultos Pequeños
A=matrix(c(0,40,70,60,0.25,0,0,0,0,0.4,0,0,0,0,0.25,0), nr=4, byrow = T)
nt <- matrix(c(400,100,250,75), ncol=1)
years=10
n.projections <- matrix(0, nrow = nrow(A), ncol = years + 1)
n.projections[, 1] <- nt
for (i in 1:years){
n.projections[, i + 1] <- A %*% n.projections[,i]
n.projections[2,] <- 0
}
n.projections
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
## [1,] 400 26000.0 6550 600 0 0 0 0 0 0 0
## [2,] 0 0.0 0 0 0 0 0 0 0 0 0
## [3,] 250 40.0 0 0 0 0 0 0 0 0 0
## [4,] 75 62.5 10 0 0 0 0 0 0 0 0
## Esta es la gráfica de la matrix si en todas las generaciones eliminamos los adultos pequenos de la poblacion
matplot(0:years, t(n.projections), type = "l", lty = 1:4,col = 1, ylab = "Stage Abundance", xlab = "Year", main = "grafica si en todas las generaciones eliminamos adulto peque?o de la poblacion")
legend("topright", legend = c("Semillas", "Adulto pequeño", "Adulto mediano", "Adulto grande"),lty = 1:4, col = 1, bty = "n")
# Para los Adultos Medianos
A=matrix(c(0,40,70,60,0.25,0,0,0,0,0.4,0,0,0,0,0.25,0), nr=4, byrow = T)
nt <- matrix(c(400,100,250,75), ncol=1)
years=10
n.projections <- matrix(0, nrow = nrow(A), ncol = years + 1)
n.projections[, 1] <- nt
for (i in 1:years){
n.projections[, i + 1] <- A %*% n.projections[,i]
n.projections[3,] <- 0
}
n.projections
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
## [1,] 400 26000.0 7750 260000.0 77500 2600000 775000 26000000 7750000
## [2,] 100 100.0 6500 1937.5 65000 19375 650000 193750 6500000
## [3,] 0 0.0 0 0.0 0 0 0 0 0
## [4,] 75 62.5 0 0.0 0 0 0 0 0
## [,10] [,11]
## [1,] 260000000 77500000
## [2,] 1937500 65000000
## [3,] 0 0
## [4,] 0 0
## Esta es la gráfica de la matrix si en todas las generaciones eliminamos los Adultos Medianos de la poblacion
matplot(0:years, t(n.projections), type = "l", lty = 1:4,col = 1, ylab = "Stage Abundance", xlab = "Year", main = "grafica si en todas las generaciones eliminamos adulto peque?o de la poblacion")
legend("topleft", legend = c("Semillas", "Adulto pequeño", "Adulto mediano", "Adulto grande"),lty = 1:4, col = 1, bty = "n")
## Puntos adicionales: Análisis de elasticidad
#El valor reproductivo es el número de individuos a nacer de una determinada edad, lo calculamos con la matriz transpuesta de la matriz de proyecciones (A)
M <- eigen(t(A))
v <- Re(M$vectors[, which.max(Re(M$values))])
v
## [1] -0.03249728 -0.45278833 -0.69324130 -0.55977063
RV <- v/v[1]
round(RV,3)
## [1] 1.000 13.933 21.332 17.225
#La importancia relativa de cada transición del ciclo de vida se hace calculando la "sensibilidad" y "elasticidad" de la matriz
#Las sensibilidades son la contribución directa de cada transición a lambda
vw.s <- v %*% t(w)
(S <- vw.s/as.numeric(v %*% w))
## [,1] [,2] [,3] [,4]
## [1,] 0.4574552 0.03283223 0.003770268 0.0002705977
## [2,] 6.3737759 0.45745524 0.052531574 0.0037702681
## [3,] 9.7585656 0.70038656 0.080428433 0.0057724666
## [4,] 7.8797360 0.56554021 0.064943440 0.0046610860
#Las elasticidades son sensibilidades ponderadas por las probabilidades de transición
elas <- (A/L1) * S
round(elas, 3)
## [,1] [,2] [,3] [,4]
## [1,] 0.000 0.377 0.076 0.005
## [2,] 0.457 0.000 0.000 0.000
## [3,] 0.000 0.080 0.000 0.000
## [4,] 0.000 0.000 0.005 0.000
# Ya que buscamos tener el maximo impacto con el mínimo esfuerzo, la transicion más adecuada para manejar esta especie invasora es el paso de semillas a adultos pequeños,
# sin embargo, cualquier manejo de las semillas implica un esfuerzo mayor y agresivo con el ambiente.
#Por tanto, es recomendable enfocarse en la transición de semilla a adulto pequeño, pero solo erradicar los adultos pequeños, ya que son estadios conspicuos más fáciles de manejar que las semillas.