stac33_a6

a. Population after 1 year

# Matrix A and initial population vector p0
A <- matrix(c(4, -2, 1, 1), nrow = 2, byrow = TRUE)
p0 <- c(30, 20)


p1 <- A %*% p0
cat("a) Population after 1 year:\n")

## a) Population after 1 year:

print(p1)

##      [,1]
## [1,]   80
## [2,]   50

After 1 year, the rabbit population increases to 80 and foxes to 50 due to reproduction and predation dynamics.

b. Population after 3 years

matrix_power <- function(M, n) {
  result <- diag(nrow(M)) 
  for (i in 1:n) {
    result <- result %*% M
  }
  return(result)
}

A3 <- matrix_power(A, 3)
p3 <- A3 %*% p0
print(p3)

##      [,1]
## [1,]  620
## [2,]  350

Over 3 years, the population continues to grow significantly.

c. Eigenvalues and eigenvectors

eig <- eigen(A)
cat("\nc) Eigenvalues:\n")

## 
## c) Eigenvalues:

print(eig$values)

## [1] 3 2

cat("Eigenvectors:\n")

## Eigenvectors:

print(eig$vectors)

##           [,1]      [,2]
## [1,] 0.8944272 0.7071068
## [2,] 0.4472136 0.7071068

The system’s dynamics are governed by eigenvalues λ = 3 and λ = 2, showing exponential growth.

d. Population after 10 years using diagonalization

c1 <- 22.36068
c2 <- 14.14214
lambda1 <- eig$values[1]
lambda2 <- eig$values[2]
v1 <- eig$vectors[, 1]
v2 <- eig$vectors[, 2]

p10 <- (lambda1^10) * c1 * v1 + (lambda2^10) * c2 * v2
cat("\nd) Population after 10 years using eigen decomposition:\n")

## 
## d) Population after 10 years using eigen decomposition:

print(round(p10))

## [1] 1191220  600730

The large values at year 10 reflect exponential growth dominated by the eigenvalue lamda = 3.

e. Create function for population at year n

pop_fn <- function(n) {
  An <- matrix_power(A, n)
  return(An %*% p0)
}

f. Use function to get population at year 3

pf3 <- pop_fn(3)
cat("\nf) Population at year 3 using function:\n")

## 
## f) Population at year 3 using function:

print(pf3)

##      [,1]
## [1,]  620
## [2,]  350

The result confirms that pop_fn() is functioning correctly. uy

g. Function for rabbits-to-foxes ratio

ratio_fn <- function(n) {
  pop <- pop_fn(n)
  if (pop[2] != 0) {
    return(pop[1] / pop[2])
  } else {
    return(NA)
  }
}
cat("\ng) Rabbit-to-fox ratio at year 3:\n")

## 
## g) Rabbit-to-fox ratio at year 3:

print(ratio_fn(3))

## [1] 1.771429

At year 3, for every fox, there are approximately 1.77 rabbits.

h. Plot rabbit-to-fox ratio over 100 years

library(ggplot2)
years <- 0:100
ratios <- sapply(years, ratio_fn)
df <- data.frame(Year = years, Ratio = ratios)

ggplot(df, aes(x = Year, y = Ratio)) +
  geom_line(color = "orange", linewidth = 1) +
  labs(title = "Rabbit-to-Fox Ratio Over Time",
       x = "Year", y = "Rabbits / Foxes") +
  theme_minimal()

The ratio converges to around 2 over time, consistent with the dominant eigenvector where rabbit:fox approximately 2:1.