if (has_sympy()) {
p <- as_sym(paste0("p", 1:3))
y <- as_sym(paste0("y", 1:3))
a <- as_sym("a")
l <- sum(y*log(p))
L <- -l + a*(sum(p) - 1)
g <- der(L, c(a, p))
sols <- solve_sys(g, c(a, p))
sol <- sols[[1L]]
sol
H <- der2(L, c(p, a))
H
H_sol <- subs_lst(H, sol)
H_sol
}
Run the code above in your browser using DataLab