Вы смотрели пакет nleqslv в R? Используя его, вы легко можете получить решение. Решение не сработает, потому что у вас есть нелинейная система уравнений, которую нельзя представить матрицами. См. Ниже описание вашей проблемы:
library(nleqslv)
eqs_fun <- function(x) {
y = numeric(2)
y[1] = pnorm(q = (x[2]+1)/5) - pnorm(q = (x[1]+1)/5) - .025
y[2] = pnorm(q = (x[2]-1)/5) - pnorm(q = (x[1]-1)/5) - .025
y
}
xstart = c(0,1)
nleqslv(xstart, eqs_fun)
Запуск этого должен дать вам следующий результат:
list(x = c(-0.159856055122849, 0.159854416799285), fvec = c(1.006390991376e-09,
-6.31374431903087e-10), termcd = 1L, message = "Function criterion near zero",
scalex = c(1, 1), nfcnt = 8L, njcnt = 1L, iter = 6L)
, и вы можете сами убедиться, что c1 = -.16, c2 = .16 будет дать вам решение!